Quick Start Guide¶
Installer¶
- Download the installer:
wget http://www.math.tau.ac.il/~haimav/libskylark-dist-0.11-Linux-x86_64.sh
- Execute it in your terminal, and follow the instructions:
bash ./libskylark-dist-0.11-Linux-x86_64.sh
Note
This should work with all fairly recent x86_64 Linux systems; please look for the string x86_64 in the output of the uname -a command.
Vagrant¶
Setting up a libSkylark environment is facilitated using Vagrant and VirtualBox.
Note
Make sure to use Vagrant >= 1.6.3! The configuration file will not be compatible with older versions.
Simply follow the instructions shown below.
Install VirtualBox (download VirtualBox)
Install Vagrant (download Vagrant).
Execute the following commands in a terminal:
git clone https://github.com/xdata-skylark/libskylark.git cd libskylark/vagrant/trusty64 vagrant up vagrant ssh
Note
The vagrant up command is downloading, building and installing libsklylark and its dependencies. This step will most likely take a while to complete.
To exit the virtual environment type exit in the shell followed by vagrant halt in order to power down the virtual machine. In order to reconnect to the virtual machine, either use the VirtualBox interface to start the virtual machine and then login with the vagrant/vagrant user/password or execute
vagrant up --no-provision
vagrant ssh
Note
Depending on your Vagrant version, omitting the --no-provision argument (after a vagrant halt) is considerably slower because the vagrant up reruns parts of the bootstrap script.
Note
The following environment variables are available after login, pointing to the libSkylark source and install directory:
SKYLARK_SRC_DIR=/home/vagrant/libskylark
SKYLARK_INSTALL_DIR=/home/vagrant/install
In the next section we provide more in depth details on how the Vagrant setup works and some trouble shooting. Additionally we provide instruction on how to install on a cluster (see Cluster of vagrant-controlled VMs) and running on AWS (see Running Vagrant on AWS).
Finally, some examples are provided in the Examples of Library Usage section.
Note
You might get a warning concerning AVX instructions when running the examples
OpenBLAS : Your OS does not support AVX instructions. OpenBLAS is using Nehalem kernels as a fallback, which may give poorer performance.
because some VM providers do not support AVX instructions yet (e.g. see: https://www.virtualbox.org/ticket/11347).
Under vagrant/trusty64 we provide an easy way to materialize an 64-bit Ubuntu-based virtual machine (release 14.04 LTS, aka Trusty Tahr) with all software dependencies required for a fresh build of libSkylark.
The Vagrantfile is the input to Vagrant, a program that facilitates the management of virtual machines. Using VirtualBox as the backend virtualization software (provider) we can launch a libSkylark environment independent of the host operating system (Windows, Mac OS X and Linux host platform require exactly the same commands).
After installing VirtualBox and Vagrant the following commands suffice
git clone https://github.com/xdata-skylark/libskylark.git
cd vagrant/trusty64
vagrant up
vagrant ssh
to get a virtual machine (VM) with libSkylark and all its dependencies. If you want to customize the installed libSkylark follow the instructions (see Building libSkylark) for building and installing libSkylark.
Note
To customize the amount of memory (4096) and CPU (75%) cap that will be used for the VM, check the Vagrant file.
The vagrant box can be stopped by issuing vagrant halt.
Note
We use clean vagrant boxes, automatically downloaded from the following urls:
Note
For a Windows XP host, PuTTy ssh client should use the default connection settings: IP= 127.0.0.1, port= 2222, username= vagrant, password= vagrant.
Hint
If you get a timeout error (waiting for the virtual machine to boot)
[...]
default: Warning: Connection timeout. Retrying...
default: Warning: Connection timeout. Retrying...
Timed out while waiting for the machine to boot. This means that
Vagrant was unable to communicate with the guest machine within
the configured ("config.vm.boot_timeout" value) time period.
during vagrant up, halt the virtual machine and vagrant up again.
Hint
In case Vagrant crashes in the bootstrap phase, the safest thing is to delete the virtual machine (halt than delete the machine in the VirtualBox GUI) and start with a fresh Vagrant build.
Command-line Usage¶
While libSkylark is designed to be a library, it also provides access to many of its feature using standalone applications. For all applications, running with the –help options reveals the command-line options.
The following example sessions assume the following commands have been performed in advance (it downloads and extracts some demo data):
wget https://www.dropbox.com/s/0ilbyeem2ihmvzw/data.tar.gz
tar -xvzf data.tar.gz
Approximate Singular Value Decomposition¶
Building libSkylark creates an executable called skylark_svd in $SKYLARK_INSTALL_DIR/bin. This executable can be used in standalone mode to compute an approximate SVD with 10 leading singular vectors as follows:
skylark_svd -k 10 --prefix usps data/usps.train
The files usps.U.txt, usps.S.txt and usps.V.txt contain the approximate SVD that was computed.
Linear Least-Squares¶
Building libSkylark creates an executable called skylark_linear in $SKYLARK_INSTALL_DIR/bin. This executable can be used in standalone mode to approximately solve linear least squares problems. Here is an example command-line:
skylark_linear data/cpu.train cpu.sol
The file cpu.sol.txt contains the computed approximate solution.
Note
Currently, only overdetermined least squares of dense matrices is supported. No regularization options are available. This will be relaxed in the future.
Learning Non-Linear Models¶
libSkylark offers two mechanisms to learn non-linear models: ADMM-based solver and kernel ridge regression solver.
Kernel Ridge Regression solver¶
In $SKYLARK_SRC_DIR/ml/skylark_krr.cpp, kernel ridge regression is setup to train and predict with a kernel based model for nonlinear classification and regression problems. The code incorporates various algorithms, with varying trade-offs in speed, memory and accuracy.
Building libSkylark creates an executable called skylark_krrl in $SKYLARK_INSTALL_DIR/bin. This executable can be used in standalone mode as follows.
- Train a classifier using the Gaussian Kernel
mpiexec -np 4 ./skylark_krr -a 1 -k 0 -g 10 -f 1000 --model model data/usps.train
- Test accuracy of the generated model and generate predictions
mpiexec -np 4 ./skylark_krr --outputfile output --model model usps.test
ADMM-based solver¶
In $SKYLARK_SRC_DIR/ml/skylark_ml.cpp, an ADMM-based solver is setup to train and predict with a randomized kernel based model for nonlinear classification and regression problems.
Building libSkylark creates an executable called skylark_ml in $SKYLARK_INSTALL_DIR/bin. This executable can be used in standalone mode as follows.
- Train an SVM with Randomized Gaussian Kernel
mpiexec -np 4 ./skylark_ml -g 10 -k 1 -l 2 -i 30 -f 1000 --trainfile data/usps.train --valfile data/usps.test --modelfile model
- Test accuracy of the generated model and generate predictions
mpiexec -np 4 ./skylark_ml --outputfile output --testfile data/usps.test --modelfile model
An output file named output.txt with the predicted labels is created.
- In the above, the entire test data is loaded to memory and the result is computed. This is the fast, but it is limited to cases where the test data is fixed. skylark_ml also supports an interactive (streaming) mode in which the software prompts for data line by line from stdin and outputs in prediction after each line is received. The mode is invoked by not passing either training data or testing data. An example for using this mode is the following:
cut -d ' ' -f 1 --complement data/usps.test | skylark_ml --modelfile model
Note
Interactive mode is much slower, since there is significant overhead that is repeated for each new line of data. In the future, batching will be supported, to avoid some of this overhead.
- It is also possible to load the model in Python and do predictions. Here is an example:
import skylark.io, skylark.ml.modeling
model = skylark.ml.modeling.LinearizedKernelModel('model')
testfile = skylark.io.libsvm('data/usps.test')
Xt, Yt = testfile.read(model.get_input_dimension())
Yp = model.predict(Xt)
Community Detection Using Seed Node¶
Building libSkylark creates an executable called skylark_community in $SKYLARK_INSTALL_DIR/bin. This executable can be used in standalone mode as follows. (Note that an interactive mode is also present.)
1. Prepare the input. Basically, the graph is described in a text file, with each row an edge. Each row has two strings that identify the node. The identifier can be arbitrary (does not have to be numeric). For example:
cat data/two_triangles
A1 A2
A2 A3
A1 A3
B1 B2
B2 B3
B1 B3
A1 B1
- Detect the community with A1 as the seed:
skylark_community --graphfile data/two_triangles --seed A1
3. The community detection algorithm is local, so it operates only on a few nodes. So, for large graphs the majority of the time for detecting a single community will be spent on just loading the file. To detect multiple communities with the graph already loaded in memory, an interactive quiet mode is provided. Here is an example:
cat data/seeds | skylark_community --graphfile data/two_triangles -i -q
Approximate Adjacency Spectral Embedding¶
Building libSkylark creates an executable called skylark_graph_se in $SKYLARK_INSTALL_DIR/bin. This executable can be used in standalone mode as follows.
1. Prepare the input. Basically, the graph is described in a text file, with each row an edge. Each row has two strings that identify the node. The identifier can be arbitrary (does not have to be numeric). For example:
cat data/two_triangles
A1 A2
A2 A3
A1 A3
B1 B2
B2 B3
B1 B3
A1 B1
- Compute a two dimensional embedding.
skylark_graph_se -k 2 data/two_triangles
3. Two files are created: out.index.txt and out.vec.txt. The vec file has the embedding vectors: each row corresponds to an index. The index file maps row index to vertices. The prefix can be called using the prefix command line option (default is out).
Cluster of vagrant-controlled VMs¶
Here is a simple approach for building a cluster of vagrant-controlled VMs; it works in a local setting but since it uses the bridged mode it should work in a real cluster environment as well. Please refer to this blog entry for more information on Networking in VirtualBox.
- In Vagrantfile set:
config.vm.network :public_network
During vagrant up choose the ethernet interface to use.
- After vagrant ssh do ifconfig, get the interface name - say
eth1 - and then for the VMs at nodes 1, 2,... do:
sudo /sbin/ifconfig eth1:vm 192.168.100.1
sudo /sbin/ifconfig eth1:vm 192.168.100.2
...
- 192.168.100.xxx will be the names to put in the “hosts” file for MPI
daemons to use; as previously noted: vagrant/vagrant is the default user/password combination for ssh.
Running Vagrant on AWS¶
This follows the instructions found at:
- https://github.com/mitchellh/vagrant-aws
- http://www.devopsdiary.com/blog/2013/05/07/automated-deployment-of-aws-ec2-instances-with-vagrant-and-puppet/
Using the provided bootstrap.sh file in combination with the aws plugin (see link above) one can deploy on AWS. First adapt the Vagrant file as shown below:
Vagrant.configure("2") do |config|
config.vm.box = "skylark"
config.vm.provision :shell, :path => "bootstrap.sh"
config.vm.provider :aws do |aws, override|
aws.access_key_id = "XYZ"
aws.secret_access_key = "AAA"
aws.keypair_name = "keynam"
aws.security_groups = ["ssh-rule"]
aws.ami = "ami-fa9cf1ca"
aws.region = "us-west-2"
aws.instance_type = "t1.micro"
override.ssh.username = "ubuntu"
override.ssh.private_key_path = "aws.pem"
end
end
Then execute:
- vagrant plugin install vagrant-aws
- vagrant box add skylark https://github.com/mitchellh/vagrant-aws/raw/master/dummy.box
- vagrant up --provider=aws
- vagrant ssh
It will take a while to compile and install everything specified in the bootstrap.sh script.
Running MPI on Hadoop/Yarn¶
libSkylark can be run on on a Yarn scheduled cluster using the blueprint outlined below.
Note
The installation requires admin privileges. Also make sure to satisfy the recommended dependencies.
1. Building MPICH2¶
Download the MPICH version specified in the dependencies (MPICH-3.1.2).
sh autogen.sh
./configure --prefix=$HOME/mpi/bin
make
make install
sudo ln -s /home/ubuntu/mpi/bin/bin/* /usr/local/bin
and copy MPI binaries to all nodes, i.e.
scp /home/ubuntu/mpi/bin/bin/* root@hadoopX:/usr/local/bin
Make sure that all the MPI executables are accessible on all the nodes.
2. Building mpich2-yarn¶
Finally we have all the dependencies ready and can continue to mpich2-yarn https://github.com/alibaba/mpich2-yarn:
cd ~/mpi
git clone https://github.com/alibaba/mpich2-yarn
cd mpich2-yarn
mvn clean package -Dmaven.test.skip=true
The configuration for Hadoop/Yarn can be found on on the Github project page. Note that the HDFS permission must be set correctly for the temporary folders used to distribute MPI executables.
Finally, an MPI job can be submitted, i.e. with
sudo -u yarn hadoop jar target/mpich2-yarn-1.0-SNAPSHOT.jar -a hellow -M 1024 -m 1024 -n 1
where we chose to run as the yarn user (double check permission for any user that runs mpich2-yarn). More examples can be found on the mpich2-yarn Github page.
Examples of Library Usage¶
Here, we provide a flavor of the library and its usage. We assume that the reader is familiar with sketching and its applications in randomized Numerical Linear algebra and Machine Learning. If not, we refer the reader to subsequent sections which provide the necessary background and references.
When libSkylark is built, executables instantiating these examples can be found under $SKYLARK_INSTALL_DIR/bin/skylark_examples and $SKYLARK_INSTALL_DIR/bin.
In addition Vagrant provisions a IPyhton notebook <http://ipython.org/notebook.html> server with access to the libSkylark Python bindings. To play with libSkylark just open your browser and direct it to localhost:6868.
Sketching¶
In $SKYLARK_SRC_DIR/examples/elemental.cpp an example is provided that illustrates sketching. Below, three types of sketches are illustrated in the highlighted lines:
- a 1D row-distributed Elemental dense matrix is sketched using the JLT (Johnson-Lindenstrauss) Transform into a local matrix
- a row-distributed Elemental dense matrix is sketched using FJLT (Fast Johnson-Lindenstrauss transform) involving faster FFT-like operations.
- a 2D block-cyclic distributed Combinatorial BLAS sparse matrix is sketched to a dense local Elemental matrix using Clarkson-Woodruff transform which involves hashing the rows.
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 132 133 134 135 136 137 138 139 140 141 142 143 144 145 146 147 148 149 150 | #include <iostream>
#include <functional>
#include <cstring>
#include <vector>
#include <utility>
#include <unordered_map>
#include <El.hpp>
#include <boost/mpi.hpp>
#define SKYLARK_NO_ANY
#include <skylark.hpp>
#include "utilities.hpp"
#include "parser.hpp"
/*******************************************/
namespace bmpi = boost::mpi;
namespace skys = skylark::sketch;
namespace skyb = skylark::base;
/*******************************************/
/* These were declared as extern in utilities.hpp --- defining it here */
int int_params [NUM_INT_PARAMETERS];
char* chr_params[NUM_CHR_PARAMETERS];
/** Typedef DistMatrix and Matrix */
typedef std::vector<int> IntContainer;
typedef std::vector<double> DblContainer;
typedef El::DistMatrix<double, El::CIRC, El::CIRC> MatrixType;
typedef El::DistMatrix<double, El::VC, El::STAR> DistMatrixType;
int main (int argc, char** argv) {
/* Initialize Elemental */
El::Initialize (argc, argv);
/* MPI sends argc and argv everywhere --- parse everywhere */
parse_parameters (argc,argv);
/* Initialize skylark */
skyb::context_t context (int_params[RAND_SEED_INDEX]);
/* Create matrices A and B */
bmpi::communicator world;
MPI_Comm mpi_world(world);
El::Grid grid (mpi_world);
El::DistMatrix<double, El::VR, El::STAR> A(grid);
El::DistMatrix<double, El::VR, El::STAR> B(grid);
/** Only randomization is supported for now */
if (0==int_params[USE_RANDOM_INDEX]) {
/** TODO: Read the entries! */
std::cout << "We don't support reading --- yet --" << std::endl;
} else {
El::Uniform (A, int_params[M_INDEX], int_params[N_INDEX]);
El::Uniform (B, int_params[M_INDEX], int_params[N_RHS_INDEX]);
}
/**
* Depending on which sketch is requested, do the sketching.
*/
if (0==strcmp("JLT", chr_params[TRANSFORM_INDEX]) ) {
if (SKETCH_LEFT == int_params[SKETCH_DIRECTION_INDEX]) {
/* 1. Create the sketching matrix */
skys::JLT_t<DistMatrixType, MatrixType> JLT (int_params[M_INDEX],
int_params[S_INDEX], context);
/* 2. Create space for the sketched matrix */
MatrixType sketch_A(int_params[S_INDEX], int_params[N_INDEX]);
/* 3. Apply the transform */
try {
JLT.apply (A, sketch_A, skys::columnwise_tag());
} catch (skylark::base::skylark_exception ex) {
SKYLARK_PRINT_EXCEPTION_DETAILS(ex);
SKYLARK_PRINT_EXCEPTION_TRACE(ex);
errno = *(boost::get_error_info<skylark::base::error_code>(ex));
std::cout << "Caught exception, exiting with error " << errno << std::endl;
std::cout << skylark_strerror(errno) << std::endl;
return errno;
}
/* 4. Print and see the result (if small enough) */
if (int_params[S_INDEX] * int_params[N_INDEX] < 100 &&
world.rank() == 0)
El::Display(sketch_A);
/** TODO: Do that same to B, and solve the system! */
} else {
std::cout << "We only have left sketching. Please retry" << std::endl;
}
} else if (0==strcmp("FJLT", chr_params[TRANSFORM_INDEX]) ) {
if (SKETCH_LEFT == int_params[SKETCH_DIRECTION_INDEX]) {
/* 1. Create the sketching matrix */
skys::FJLT_t<DistMatrixType, MatrixType> FJLT (int_params[M_INDEX],
int_params[S_INDEX], context);
/* 2. Create space for the sketched matrix */
MatrixType sketch_A(int_params[S_INDEX], int_params[N_INDEX]);
/* 3. Apply the transform */
FJLT.apply (A, sketch_A, skys::columnwise_tag());
/* 4. Print and see the result */
if (int_params[S_INDEX] * int_params[M_INDEX] < 100 &&
world.rank() == 0)
El::Display(sketch_A);
/** TODO: Do that same to B, and solve the system! */
}
} else if (0==strcmp("CWT", chr_params[TRANSFORM_INDEX]) ) {
if (SKETCH_LEFT == int_params[SKETCH_DIRECTION_INDEX]) {
/* 1. Create the sketching matrix */
skys::CWT_t<DistMatrixType, MatrixType>
Sparse (int_params[M_INDEX], int_params[S_INDEX], context);
/* 2. Create space for the sketched matrix */
MatrixType sketch_A(int_params[S_INDEX], int_params[N_INDEX]);
/* 3. Apply the transform */
SKYLARK_BEGIN_TRY()
Sparse.apply (A, sketch_A, skys::columnwise_tag());
SKYLARK_END_TRY()
SKYLARK_CATCH_AND_RETURN_ERROR_CODE();
/* 4. Print and see the result */
if (int_params[S_INDEX] * int_params[M_INDEX] < 100 &&
world.rank() == 0)
El::Display(sketch_A);
/** TODO: Do that same to B, and solve the system! */
}
} else {
std::cout << "We only have JLT/FJLT/Sparse sketching. Please retry" <<
std::endl;
}
El::Finalize();
return 0;
}
|
For sketching Elemental and CombBLAS matrices via Python bindings, run, for example, the following:
mpiexec -np 4 python $SKYLARK_SRC_DIR/python-skylark/skylark/examples/example_sketch.py
mpiexec -np 4 python $SKYLARK_SRC_DIR/python-skylark/skylark/examples/example_sparse_sketch.py
Note
Currently, we have stable python bindings only for sketching layer. This too might change with newer versions of Elemental.
Fast Least Square Regression¶
In $SKYLARK_SRC_DIR/examples/least_squares.cpp examples are provided on how least squares problems can be solved faster using sketching. One approach is of the flavor of sketch-and-solve geared towards lower-precision solutions, while the other approach uses sketching to construct a preconditioner to obtain high-precision solutions faster.
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 132 133 | #include <iostream>
#include <El.hpp>
#include <boost/mpi.hpp>
#include <boost/format.hpp>
#define SKYLARK_NO_ANY
#include <skylark.hpp>
const int m = 50000;
const int n = 500;
template<typename MatrixType, typename RhsType, typename SolType>
void check_solution(const MatrixType &A, const RhsType &b, const SolType &x,
const RhsType &r0,
double &res, double &resAtr, double &resFac) {
RhsType r(b);
skylark::base::Gemv(El::NORMAL, -1.0, A, x, 1.0, r);
res = skylark::base::Nrm2(r);
SolType Atr(x.Height(), x.Width(), x.Grid());
skylark::base::Gemv(El::TRANSPOSE, 1.0, A, r, 0.0, Atr);
resAtr = skylark::base::Nrm2(Atr);
skylark::base::Axpy(-1.0, r0, r);
RhsType dr(b);
skylark::base::Axpy(-1.0, r0, dr);
resFac = skylark::base::Nrm2(r) / skylark::base::Nrm2(dr);
}
template<typename MatrixType, typename RhsType, typename SolType>
void experiment() {
typedef MatrixType matrix_type;
typedef RhsType rhs_type;
typedef SolType sol_type;
double res, resAtr, resFac;
boost::mpi::communicator world;
int rank = world.rank();
skylark::base::context_t context(23234);
// Setup problem and righthand side
// Using Skylark's uniform generator (as opposed to Elemental's)
// will insure the same A and b are generated regardless of the number
// of processors.
matrix_type A, b;
skylark::base::UniformMatrix(A, m, n, context);
skylark::base::UniformMatrix(b, m, 1, context);
sol_type x(n,1);
rhs_type r(b);
boost::mpi::timer timer;
double telp;
// Solve using Elemental. Note: Elemental only supports [MC,MR]...
El::DistMatrix<double> A1 = A, b1 = b, x1;
timer.restart();
El::LeastSquares(El::NORMAL, A1, b1, x1);
telp = timer.elapsed();
x = x1;
check_solution(A, b, x, r, res, resAtr, resFac);
if (rank == 0)
std::cout << "Elemental:\t\t\t||r||_2 = "
<< boost::format("%.2f") % res
<< "\t\t\t\t\t\t\t||A' * r||_2 = " << boost::format("%.2e") % resAtr
<< "\t\tTime: " << boost::format("%.2e") % telp << " sec"
<< std::endl;
double res_opt = res;
// The following computes the optimal residual (r^\star in the logs)
skylark::base::Gemv(El::NORMAL, -1.0, A, x, 1.0, r);
#if SKYLARK_HAVE_FFTW || SKYLARK_HAVE_FFTWF
// Solve using Sylark
timer.restart();
skylark::nla::FasterLeastSquares(El::NORMAL, A, b, x, context);
telp = timer.elapsed();
check_solution(A, b, x, r, res, resAtr, resFac);
if (rank == 0)
std::cout << "Skylark:\t\t\t||r||_2 = "
<< boost::format("%.2f") % res
<< " (x " << boost::format("%.5f") % (res / res_opt) << ")"
<< "\t||r - r*||_2 / ||b - r*||_2 = " << boost::format("%.2e") % resFac
<< "\t||A' * r||_2 = " << boost::format("%.2e") % resAtr
<< "\t\tTime: " << boost::format("%.2e") % telp << " sec"
<< std::endl;
// Approximately solve using Sylark
timer.restart();
skylark::nla::ApproximateLeastSquares(El::NORMAL, A, b, x, context);
telp = timer.elapsed();
check_solution(A, b, x, r, res, resAtr, resFac);
if (rank == 0)
std::cout << "Skylark (approximate):\t\t||r||_2 = "
<< boost::format("%.2f") % res
<< " (x " << boost::format("%.5f") % (res / res_opt) << ")"
<< "\t||r - r*||_2 / ||b - r*||_2 = " << boost::format("%.2e") % resFac
<< "\t||A' * r||_2 = " << boost::format("%.2e") % resAtr
<< "\t\tTime: " << boost::format("%.2e") % telp << " sec"
<< std::endl;
#else
std::cout << "You need to have Skylark supporting FFTW or FFTWF "
<< "to solve with skylark least_squares.cpp"
<< std::endl;
#endif
}
int main(int argc, char** argv) {
El::Initialize(argc, argv);
boost::mpi::communicator world;
int rank = world.rank();
if (rank == 0)
std::cout << "Matrix: [VC,STAR], Rhs: [VC,STAR], Sol: [STAR,STAR]\n\n";
experiment<El::DistMatrix<double, El::VC, El::STAR>,
El::DistMatrix<double, El::VC, El::STAR>,
El::DistMatrix<double, El::STAR, El::STAR> > ();
if (rank == 0)
std::cout << "\nMatrix: [MC,MR], Rhs: [MC,MR], Sol: [MC,MR]\n\n";
experiment<El::DistMatrix<double>, El::DistMatrix<double>,
El::DistMatrix<double> >();
return 0;
}
|
For an example of sketch-and-solve regression in Python, run:
wget http://vikas.sindhwani.org/data.tar.gz
tar -xvzf data.tar.gz
mpiexec -np 4 python $SKYLARK_SRC_DIR/python-skylark/skylark/examples/sketch_solve.py data/usps.train data/usps.test
SVD¶
In $SKYLARK_SRC_DIR/nla/skylark_svd.cpp an example is provided that illustrates randomized singular value decompositions.
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 132 133 134 135 136 137 138 139 140 141 142 143 144 145 146 147 148 149 150 151 152 153 154 155 156 157 158 159 160 161 162 163 164 165 166 167 168 169 170 171 172 173 174 175 176 177 178 179 180 181 182 183 184 185 186 187 188 189 190 191 192 193 194 195 196 197 198 199 200 201 202 203 204 205 206 207 208 209 210 211 212 213 214 215 216 217 218 219 220 221 222 223 224 225 226 227 228 229 230 231 232 233 234 235 236 237 238 239 240 241 242 243 244 245 246 247 248 249 250 251 252 253 254 255 256 257 258 259 260 261 262 263 264 265 266 267 268 269 270 271 272 273 274 275 276 277 278 279 280 281 282 283 284 285 286 287 288 289 290 291 292 293 294 295 296 297 298 299 300 301 302 303 304 305 306 307 308 309 310 311 312 313 314 315 316 317 318 319 320 321 322 323 324 325 326 327 328 329 330 331 332 333 334 335 336 337 338 339 340 341 342 343 344 345 346 347 348 349 350 351 352 353 354 355 356 357 358 359 360 361 362 363 364 365 366 367 368 369 370 371 372 373 374 375 376 377 378 379 380 381 382 383 384 385 386 387 388 389 390 391 392 393 394 395 396 397 398 399 400 401 402 403 404 405 406 407 408 409 410 411 412 413 414 415 416 417 418 419 420 421 422 423 424 425 426 427 428 429 430 431 432 433 434 435 436 437 438 439 440 441 442 443 444 445 446 447 448 449 450 451 452 453 454 455 456 457 458 459 460 461 462 463 464 465 466 467 468 469 470 471 472 473 474 475 476 477 | #include <El.hpp>
#include <boost/program_options.hpp>
#include <boost/mpi.hpp>
#include <boost/format.hpp>
#define SKYLARK_NO_ANY
#include <skylark.hpp>
#include <iostream>
namespace bpo = boost::program_options;
enum file_types {
LIBSVM,
ARC_LIST
};
template<typename InputType, typename FactorType, typename UType = FactorType,
typename YType = FactorType>
void execute(bool directory, const std::string &fname,
const std::string &hdfs, int port, const std::vector<int> &profile, int k,
const skylark::nla::approximate_svd_params_t ¶ms,
const std::string &prefix,
skylark::base::context_t &context) {
boost::mpi::communicator world;
int rank = world.rank();
InputType A;
FactorType S, V;
UType U;
YType Y;
boost::mpi::timer timer;
if (profile.empty()) {
// Load A and Y (Y is thrown away)
if (rank == 0) {
std::cout << "Reading the matrix... ";
std::cout.flush();
timer.restart();
}
if (!hdfs.empty()) {
# if SKYLARK_HAVE_LIBHDFS
hdfsFS fs;
if (rank == 0)
fs = hdfsConnect(hdfs.c_str(), port);
if (directory)
SKYLARK_THROW_EXCEPTION(skylark::base::io_exception() <<
skylark::base::error_msg(
"HDFS directory reading not yet supported."))
else
skylark::utility::io::ReadLIBSVM(
fs, fname, A, Y, skylark::base::ROWS);
# else
SKYLARK_THROW_EXCEPTION(skylark::base::io_exception() <<
skylark::base::error_msg("Install libhdfs for HDFS support!"));
# endif
} else {
if (directory)
skylark::utility::io::ReadDirLIBSVM(
fname, A, Y, skylark::base::ROWS);
else
skylark::utility::io::ReadLIBSVM(fname, A, Y, skylark::base::ROWS);
}
Y.Empty();
} else {
if (rank == 0) {
std::cout << "Generating random matrix... ";
std::cout.flush();
timer.restart();
}
skylark::base::UniformMatrix(A, profile[0], profile[1], context);
}
if (rank == 0)
std::cout <<"took " << boost::format("%.2e") % timer.elapsed()
<< " sec\n";
/* Compute approximate SVD */
if (rank == 0) {
std::cout << "Computing approximate SVD...";
std::cout.flush();
timer.restart();
}
skylark::nla::ApproximateSVD(A, U, S, V, k, context, params);
if (rank == 0)
std::cout <<"Took " << boost::format("%.2e") % timer.elapsed()
<< " sec\n";
/* Write results */
if (rank == 0) {
std::cout << "Writing results...";
std::cout.flush();
timer.restart();
}
El::Write(U, prefix + ".U", El::ASCII);
El::Write(S, prefix + ".S", El::ASCII);
El::Write(V, prefix + ".V", El::ASCII);
if (rank == 0)
std::cout <<"took " << boost::format("%.2e") % timer.elapsed()
<< " sec\n";
}
template<typename InputType, typename FactorType, typename UType = FactorType,
typename YType = FactorType>
void execute_sym(bool directory, const std::string &fname,
const std::string& ftype,
const std::string &hdfs, int port,
const std::vector<int> &profile, bool lower, int k,
const skylark::nla::approximate_svd_params_t ¶ms,
const std::string &prefix,
skylark::base::context_t &context) {
boost::mpi::communicator world;
int rank = world.rank();
InputType A;
FactorType S, V;
YType Y;
boost::mpi::timer timer;
if (profile.empty()) {
// Load A and Y (Y is thrown away)
if (rank == 0) {
std::cout << "Reading the matrix... ";
std::cout.flush();
timer.restart();
}
if (!hdfs.empty()) {
# if SKYLARK_HAVE_LIBHDFS
hdfsFS fs;
if (rank == 0)
fs = hdfsConnect(hdfs.c_str(), port);
if (directory)
SKYLARK_THROW_EXCEPTION(skylark::base::io_exception() <<
skylark::base::error_msg(
"HDFS directory reading not yet supported."))
else
skylark::utility::io::ReadLIBSVM(
fs, fname, A, Y, skylark::base::ROWS);
# else
SKYLARK_THROW_EXCEPTION(skylark::base::io_exception() <<
skylark::base::error_msg("Install libhdfs for HDFS support!"));
# endif
} else {
// FIXME: ugly, fix options
if (ftype.compare("ARC_LIST") == 0) {
skylark::utility::io::ReadArcList(fname, A, world, true);
} else {
if (directory)
skylark::utility::io::ReadDirLIBSVM(
fname, A, Y, skylark::base::ROWS);
else
skylark::utility::io::ReadLIBSVM(
fname, A, Y, skylark::base::ROWS);
}
}
Y.Empty();
} else {
SKYLARK_THROW_EXCEPTION(skylark::base::unsupported_base_operation() <<
skylark::base::error_msg("Uniform symmetric matrix generating not supported yet."));
}
if (rank == 0)
std::cout << "took " << boost::format("%.2e") % timer.elapsed()
<< " sec\n";
/* Compute approximate SVD */
if (rank == 0) {
std::cout << "Computing approximate SVD...";
std::cout.flush();
timer.restart();
}
skylark::nla::ApproximateSymmetricSVD(lower ? El::LOWER : El::UPPER,
A, V, S, k, context, params);
if (rank == 0)
std::cout <<"Took " << boost::format("%.2e") % timer.elapsed()
<< " sec\n";
/* Write results */
if (rank == 0) {
std::cout << "Writing results...";
std::cout.flush();
timer.restart();
}
El::Write(S, prefix + ".S", El::ASCII);
El::Write(V, prefix + ".V", El::ASCII);
if (rank == 0)
std::cout <<"took " << boost::format("%.2e") % timer.elapsed()
<< " sec\n";
world.barrier();
}
int main(int argc, char* argv[]) {
El::Initialize(argc, argv);
boost::mpi::communicator world;
int rank = world.rank();
int size = world.size();
int seed, k, powerits, port;
std::string fname, ftype, prefix, hdfs;
bool as_symmetric, as_sparse, skipqr, use_single, lower, directory;
int oversampling_ratio, oversampling_additive;
std::vector<int> profile;
// Parse options
bpo::options_description desc("Options");
desc.add_options()
("help,h", "produce a help message")
("inputfile",
bpo::value<std::string>(&fname),
"Input file to run approximate SVD on (default libsvm format).")
("filetype",
bpo::value<std::string>(&ftype),
"Input file type (LIBSVM or ARC_LIST).")
("directory,d", "Whether inputfile is a directory of files whose"
" concatination is the input.")
("seed,s",
bpo::value<int>(&seed)->default_value(38734),
"Seed for random number generation. OPTIONAL.")
("hdfs",
bpo::value<std::string>(&hdfs)->default_value(""),
"If not empty, will assume file is in an HDFS. "
"Parameter is filesystem name.")
("port",
bpo::value<int>(&port)->default_value(0),
"For HDFS: port to use.")
("rank,k",
bpo::value<int>(&k)->default_value(6),
"Target rank. OPTIONAL.")
("powerits,i",
bpo::value<int>(&powerits)->default_value(2),
"Number of power iterations. OPTIONAL.")
("skipqr", "Whether to skip QR in each iteration. Higher than one power"
" iterations is not recommended in this mode.")
("ratio,r",
bpo::value<int>(&oversampling_ratio)->default_value(2),
"Ratio of oversampling of rank. OPTIONAL.")
("additive,a",
bpo::value<int>(&oversampling_additive)->default_value(0),
"Additive factor for oversampling of rank. OPTIONAL.")
("symmetric", "Whether to treat the matrix as symmetric. "
"Only upper part will be accessed, unless --lower is used.")
("lower", "For symmetric matrix, access only lower part (upper is default).")
("sparse", "Whether to load the matrix as a sparse one.")
("single", "Whether to use single precision instead of double.")
("profile",
bpo::value<std::vector<int> >()->multitoken(),
"Generate random matrix and run on it (for profiling)."
"Requires specification of height and width. OPTIONAL.")
("prefix",
bpo::value<std::string>(&prefix)->default_value("out"),
"Prefix for output files (prefix.U.txt, prefix.S.txt"
" and prefix.V.txt. OPTIONAL.");
bpo::positional_options_description positional;
positional.add("inputfile", 1);
bpo::variables_map vm;
try {
bpo::store(bpo::command_line_parser(argc, argv)
.options(desc).positional(positional).run(), vm);
if (vm.count("help")) {
if (rank == 0) {
std::cout << "Usage: " << argv[0]
<< " [options] input-file-name" << std::endl;
std::cout << desc;
}
world.barrier();
return 0;
}
if (vm["profile"].empty() && !vm.count("inputfile")) {
if (rank == 0)
std::cout << "Input file is required." << std::endl;
world.barrier();
return -1;
}
if (!vm["profile"].empty()) {
if (vm["profile"].as<std::vector<int> >().size() < 2) {
if (rank == 0)
std::cout << "Please specify height and width for --profile."
<< std::endl;
world.barrier();
return -1;
} else
profile = vm["profile"].as<std::vector<int> >();
}
bpo::notify(vm);
as_symmetric = vm.count("symmetric");
as_sparse = vm.count("sparse");
skipqr = vm.count("skipqr");
use_single = vm.count("single");
directory = vm.count("directory");
lower = vm.count("lower");
} catch(bpo::error& e) {
if (rank == 0) {
std::cerr << e.what() << std::endl;
std::cerr << desc << std::endl;
}
world.barrier();
return -1;
}
skylark::base::context_t context(seed);
skylark::nla::approximate_svd_params_t params;
params.skip_qr = skipqr;
params.num_iterations = powerits;
params.oversampling_ratio = oversampling_ratio;
params.oversampling_additive = oversampling_additive;
SKYLARK_BEGIN_TRY()
if (size == 1) {
if (!as_symmetric) {
if (use_single) {
if (as_sparse)
execute<skylark::base::sparse_matrix_t<float>,
El::Matrix<float> >(
directory, fname, hdfs,
port, profile, k, params, prefix, context);
else
execute<El::Matrix<float>,
El::Matrix<float> >(directory, fname, hdfs,
port, profile, k, params, prefix, context);
} else {
if (as_sparse)
execute<skylark::base::sparse_matrix_t<double>,
El::Matrix<double> >(
directory, fname, hdfs,
port, profile, k, params, prefix, context);
else
execute<El::Matrix<double>,
El::Matrix<double> >(directory, fname, hdfs,
port, profile, k, params, prefix, context);
}
} else {
if (use_single) {
if (as_sparse)
execute_sym<skylark::base::sparse_matrix_t<float>,
El::Matrix<float> >(
directory, fname, ftype,
hdfs, port, profile, lower, k, params, prefix,
context);
else
execute_sym<El::Matrix<float>,
El::Matrix<float> >(directory, fname, ftype,
hdfs, port, profile, lower, k, params, prefix,
context);
} else {
if (as_sparse)
execute_sym<skylark::base::sparse_matrix_t<double>,
El::Matrix<double>,
El::Matrix<double> >(
directory, fname, ftype, hdfs, port,
profile, lower, k, params, prefix, context);
else
execute_sym<El::Matrix<double>,
El::Matrix<double> >(directory, fname,
ftype, hdfs, port, profile, lower, k, params,
prefix, context);
}
}
} else {
if (!as_symmetric) {
if (use_single) {
if (as_sparse)
execute<skylark::base::sparse_vc_star_matrix_t<float>,
El::DistMatrix<float, El::STAR, El::STAR>,
El::DistMatrix<float, El::VC, El::STAR>,
El::DistMatrix<float, El::VC, El::STAR> >(
directory, fname, hdfs,
port, profile, k, params, prefix, context);
else
execute<El::DistMatrix<float>,
El::DistMatrix<float> >(directory, fname, hdfs,
port, profile, k, params, prefix, context);
} else {
if (as_sparse)
execute<skylark::base::sparse_vc_star_matrix_t<double>,
El::DistMatrix<double, El::STAR, El::STAR>,
El::DistMatrix<double, El::VC, El::STAR>,
El::DistMatrix<double, El::VC, El::STAR> >(
directory, fname, hdfs,
port, profile, k, params, prefix, context);
else
execute<El::DistMatrix<double>,
El::DistMatrix<double> >(directory, fname, hdfs,
port, profile, k, params, prefix, context);
}
} else {
if (use_single) {
if (as_sparse)
execute_sym<skylark::base::sparse_vc_star_matrix_t<float>,
El::DistMatrix<float, El::VC, El::STAR>,
El::DistMatrix<float, El::VC, El::STAR> >(
directory, fname, ftype,
hdfs, port, profile, lower, k, params, prefix,
context);
else
execute_sym<El::DistMatrix<float>,
El::DistMatrix<float> >(directory, fname, ftype,
hdfs, port, profile, lower, k, params, prefix,
context);
} else {
if (as_sparse)
execute_sym<skylark::base::sparse_vc_star_matrix_t<double>,
El::DistMatrix<double, El::VC, El::STAR>,
El::DistMatrix<double, El::VC, El::STAR> >(
directory, fname, ftype, hdfs, port,
profile, lower,
k, params, prefix, context);
else
execute_sym<El::DistMatrix<double>,
El::DistMatrix<double> >(directory, fname,
ftype, hdfs, port, profile, lower, k, params,
prefix, context);
}
}
}
SKYLARK_END_TRY() SKYLARK_CATCH_AND_PRINT((rank == 0))
El::Finalize();
return 0;
}
|