]> https://gitweb.dealii.org/ - dealii.git/commitdiff
enable chunking for save/load of ScaLAPACKMatrix 5867/head
authorBenjamin Brands <benjamin.brands@fau.de>
Tue, 6 Feb 2018 06:59:04 +0000 (07:59 +0100)
committerBenjamin Brands <benjamin.brands@fau.de>
Wed, 7 Feb 2018 18:54:26 +0000 (19:54 +0100)
include/deal.II/lac/scalapack.h
source/lac/scalapack.cc
tests/scalapack/scalapack_10_b.cc [new file with mode: 0644]
tests/scalapack/scalapack_10_b.with_hdf5=on.mpirun=1.output [new file with mode: 0644]
tests/scalapack/scalapack_10_b.with_hdf5=on.mpirun=11.output [new file with mode: 0644]
tests/scalapack/scalapack_10_b.with_hdf5=on.mpirun=4.output [new file with mode: 0644]
tests/scalapack/scalapack_10_b.with_hdf5=on.mpirun=5.output [new file with mode: 0644]

index 42468c81fc958e2c1b71e37253ee8be366564aac..9686e91918c704e77bb97a43adbf0c774f58a2c4 100644 (file)
@@ -191,6 +191,7 @@ public:
 
   /**
    * Stores the distributed matrix in @p filename using HDF5.
+   *
    * In case that deal.II was built without HDF5
    * a call to this function will cause an exception to be thrown.
    *
@@ -199,8 +200,15 @@ public:
    * internally the distributed matrix is copied to one process, which
    * does the output. Therefore, the matrix has to fit into the memory
    * of one process.
+   *
+   * To tweak the I/O performance, especially for parallel I/O, the user may define the optional parameter @p chunk_size.
+   * All MPI processes need to call the function with the same value.
+   * The matrix is written in chunks to the file, therefore the properties of the system define the optimal chunk size.
+   * Internally, HDF5 splits the matrix into <tt>chunk_size.first</tt> x <tt>chunk_size.second</tt> sized blocks,
+   * with <tt>chunk_size.first</tt> being the number of rows of a chunk and <tt>chunk_size.second</tt> the number of columns.
    */
-  void save(const char *filename) const;
+  void save(const char *filename,
+            const std::pair<unsigned int,unsigned int> &chunk_size=std::make_pair(numbers::invalid_unsigned_int,numbers::invalid_unsigned_int)) const;
 
   /**
    * Loads the distributed matrix from file @p filename using HDF5.
@@ -417,7 +425,8 @@ private:
    * Stores the distributed matrix in @p filename
    * using serial routines
    */
-  void save_serial(const char *filename) const;
+  void save_serial(const char *filename,
+                   const std::pair<unsigned int,unsigned int> &chunk_size) const;
 
   /*
    * Loads the distributed matrix from file @p filename
@@ -429,7 +438,8 @@ private:
    * Stores the distributed matrix in @p filename
    * using parallel routines
    */
-  void save_parallel(const char *filename) const;
+  void save_parallel(const char *filename,
+                     const std::pair<unsigned int,unsigned int> &chunk_size) const;
 
   /*
    * Loads the distributed matrix from file @p filename
index 55d592b31a9e505d00a7156320de694d930dcecd..f87838f1fa96cc5643a935885f2db884d79d945d 100644 (file)
@@ -914,19 +914,34 @@ NumberType ScaLAPACKMatrix<NumberType>::norm(const char type) const
 
 
 template <typename NumberType>
-void ScaLAPACKMatrix<NumberType>::save(const char *filename) const
+void ScaLAPACKMatrix<NumberType>::save(const char *filename,
+                                       const std::pair<unsigned int,unsigned int> &chunk_size) const
 {
 #ifndef DEAL_II_WITH_HDF5
   (void)filename;
+  (void)chunk_size;
   AssertThrow(false, ExcMessage ("HDF5 support is disabled."));
 #else
+
+  Assert((chunk_size.first <= (unsigned int)n_rows) && (chunk_size.first>0),ExcIndexRange(chunk_size.first,1,n_rows+1));
+  Assert((chunk_size.second <= (unsigned int)n_columns) && (chunk_size.second>0),ExcIndexRange(chunk_size.second,1,n_columns+1));
+
+  std::pair<unsigned int,unsigned int> chunks_size_ = chunk_size;
+
+  if (chunks_size_.first==numbers::invalid_unsigned_int || chunks_size_.second==numbers::invalid_unsigned_int)
+    {
+      // default: store the matrix in chunks of columns
+      chunks_size_.first = n_rows;
+      chunks_size_.second = 1;
+    }
+
 #  ifdef H5_HAVE_PARALLEL
   //implementation for configurations equipped with a parallel file system
-  save_parallel(filename);
+  save_parallel(filename,chunks_size_);
 
 #  else
   //implementation for configurations with no parallel file system
-  save_serial(filename);
+  save_serial(filename,chunks_size_);
 
 #  endif
 #endif
@@ -935,10 +950,12 @@ void ScaLAPACKMatrix<NumberType>::save(const char *filename) const
 
 
 template <typename NumberType>
-void ScaLAPACKMatrix<NumberType>::save_serial(const char *filename) const
+void ScaLAPACKMatrix<NumberType>::save_serial(const char *filename,
+                                              const std::pair<unsigned int,unsigned int> &chunk_size) const
 {
 #  ifndef DEAL_II_WITH_HDF5
   (void)filename;
+  (void)chunk_size;
   Assert(false,ExcInternalError());
 #  else
 
@@ -964,6 +981,15 @@ void ScaLAPACKMatrix<NumberType>::save_serial(const char *filename) const
       // create a new file using default properties
       hid_t file_id = H5Fcreate(filename, H5F_ACC_TRUNC, H5P_DEFAULT, H5P_DEFAULT);
 
+      // modify dataset creation properties, i.e. enable chunking
+      hsize_t chunk_dims[2];
+      //revert order of rows and columns as ScaLAPACK uses column-major ordering
+      chunk_dims[0] = chunk_size.second;
+      chunk_dims[1] = chunk_size.first;
+      hid_t property = H5Pcreate (H5P_DATASET_CREATE);
+      status = H5Pset_chunk (property, 2, chunk_dims);
+      AssertThrow(status >= 0, ExcIO());
+
       // create the data space for the dataset
       hsize_t dims[2];
       //change order of rows and columns as ScaLAPACKMatrix uses column major ordering
@@ -971,11 +997,11 @@ void ScaLAPACKMatrix<NumberType>::save_serial(const char *filename) const
       dims[1] = n_rows;
       hid_t dataspace_id = H5Screate_simple(2, dims, nullptr);
 
-      // create the dataset
+      // create the dataset within the file using chunk creation properties
       hid_t type_id = hdf5_type_id(&tmp.values[0]);
       hid_t dataset_id = H5Dcreate2(file_id, "/matrix",
                                     type_id, dataspace_id,
-                                    H5P_DEFAULT, H5P_DEFAULT, H5P_DEFAULT);
+                                    H5P_DEFAULT, property, H5P_DEFAULT);
 
       // write the dataset
       status = H5Dwrite(dataset_id, type_id,
@@ -991,6 +1017,10 @@ void ScaLAPACKMatrix<NumberType>::save_serial(const char *filename) const
       status = H5Sclose(dataspace_id);
       AssertThrow(status >= 0, ExcIO());
 
+      // release the creation property
+      status = H5Pclose (property);
+      AssertThrow(status >= 0, ExcIO());
+
       // close the file.
       status = H5Fclose(file_id);
       AssertThrow(status >= 0, ExcIO());
@@ -1001,10 +1031,12 @@ void ScaLAPACKMatrix<NumberType>::save_serial(const char *filename) const
 
 
 template <typename NumberType>
-void ScaLAPACKMatrix<NumberType>::save_parallel(const char *filename) const
+void ScaLAPACKMatrix<NumberType>::save_parallel(const char *filename,
+                                                const std::pair<unsigned int,unsigned int> &chunk_size) const
 {
 #  ifndef DEAL_II_WITH_HDF5
   (void)filename;
+  (void)chunk_size;
   Assert(false,ExcInternalError());
 #  else
 
@@ -1048,12 +1080,20 @@ void ScaLAPACKMatrix<NumberType>::save_parallel(const char *filename) const
 
   hid_t filespace = H5Screate_simple(2, dims, nullptr);
 
-  // create the dataset with default properties and close filespace
+  // create the chunked dataset with default properties and close filespace
+  hsize_t chunk_dims[2];
+  //revert order of rows and columns as ScaLAPACK uses column-major ordering
+  chunk_dims[0] = chunk_size.second;
+  chunk_dims[1] = chunk_size.first;
+  plist_id = H5Pcreate(H5P_DATASET_CREATE);
+  H5Pset_chunk(plist_id, 2, chunk_dims);
   hid_t type_id = hdf5_type_id(data);
   hid_t dset_id = H5Dcreate2(file_id, "/matrix", type_id,
-                             filespace, H5P_DEFAULT, H5P_DEFAULT, H5P_DEFAULT);
+                             filespace, H5P_DEFAULT, plist_id, H5P_DEFAULT);
   status = H5Sclose(filespace);
   AssertThrow(status >= 0, ExcIO());
+  status = H5Pclose(plist_id);
+  AssertThrow(status >= 0, ExcIO());
 
   // gather the number of local rows and columns from all processes
   std::vector<int> proc_n_local_rows(n_mpi_processes), proc_n_local_columns(n_mpi_processes);
diff --git a/tests/scalapack/scalapack_10_b.cc b/tests/scalapack/scalapack_10_b.cc
new file mode 100644 (file)
index 0000000..3e09b51
--- /dev/null
@@ -0,0 +1,96 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2017 by the deal.II authors
+//
+// This file is part of the deal.II library.
+//
+// The deal.II library is free software; you can use it, redistribute
+// it, and/or modify it under the terms of the GNU Lesser General
+// Public License as published by the Free Software Foundation; either
+// version 2.1 of the License, or (at your option) any later version.
+// The full text of the license can be found in the file LICENSE at
+// the top level of the deal.II distribution.
+//
+// ---------------------------------------------------------------------
+
+#include "../tests.h"
+#include "../lapack/create_matrix.h"
+
+// test serial saving and loading of distributed ScaLAPACKMatrices with prescribed chunk sizes
+
+#include <deal.II/base/logstream.h>
+#include <deal.II/base/utilities.h>
+#include <deal.II/base/conditional_ostream.h>
+#include <deal.II/base/timer.h>
+#include <deal.II/base/multithread_info.h>
+
+#include <deal.II/lac/scalapack.h>
+
+#include <fstream>
+#include <iostream>
+#include <cstdio>
+
+
+template <typename NumberType>
+void test(const std::pair<unsigned int,unsigned int> &size, const unsigned int block_size, const std::pair<unsigned int,unsigned int> &chunk_size)
+{
+  const std::string filename ("scalapck_10_b_test.h5");
+
+  MPI_Comm mpi_communicator(MPI_COMM_WORLD);
+  const unsigned int this_mpi_process(Utilities::MPI::this_mpi_process(mpi_communicator));
+  ConditionalOStream pcout (std::cout, (this_mpi_process ==0));
+
+  FullMatrix<NumberType> full(size.first,size.second);
+  create_random(full);
+
+  //create 2d process grid
+  std::shared_ptr<Utilities::MPI::ProcessGrid> grid = std::make_shared<Utilities::MPI::ProcessGrid>(mpi_communicator,size.first,
+                                                      size.second,block_size,block_size);
+
+  ScaLAPACKMatrix<NumberType> scalapack_matrix(size.first,size.second,grid,block_size,block_size);
+  ScaLAPACKMatrix<NumberType> scalapack_matrix_copy(size.first,size.second,grid,block_size,block_size);
+
+  scalapack_matrix = full;
+  scalapack_matrix.save(filename.c_str(),chunk_size);
+  scalapack_matrix_copy.load(filename.c_str());
+
+  FullMatrix<NumberType> copy(size.first,size.second);
+  scalapack_matrix_copy.copy_to(copy);
+  copy.add(-1,full);
+
+  pcout << size.first << "x" << size.second << " & "
+        << block_size << " & "
+        << chunk_size.first << "x" << chunk_size.second << std::endl;
+  AssertThrow(copy.frobenius_norm()<1e-12,ExcInternalError());
+  std::remove(filename.c_str());
+}
+
+
+
+int main (int argc,char **argv)
+{
+  Utilities::MPI::MPI_InitFinalize mpi_initialization(argc, argv, numbers::invalid_unsigned_int);
+
+  std::vector<std::pair<unsigned int,unsigned int>> sizes;
+  sizes.push_back(std::make_pair(100,75));
+  sizes.push_back(std::make_pair(200,225));
+  sizes.push_back(std::make_pair(300,250));
+
+  const std::vector<unsigned int> block_sizes = {{1,16,32}};
+
+  std::vector<std::pair<unsigned int,unsigned int>> chunk_sizes;
+  chunk_sizes.push_back(std::make_pair(1,1));
+  chunk_sizes.push_back(std::make_pair(10,10));
+  chunk_sizes.push_back(std::make_pair(50,50));
+  chunk_sizes.push_back(std::make_pair(100,75));
+
+  for (unsigned int i=0; i<sizes.size(); ++i)
+    for (unsigned int j=0; j<block_sizes.size(); ++j)
+      for (unsigned int k=0; k<chunk_sizes.size(); ++k)
+        test<double>(sizes[i],block_sizes[j],chunk_sizes[k]);
+
+  for (unsigned int i=0; i<sizes.size(); ++i)
+    for (unsigned int j=0; j<block_sizes.size(); ++j)
+      for (unsigned int k=0; k<chunk_sizes.size(); ++k)
+        test<float>(sizes[i],block_sizes[j],chunk_sizes[k]);
+}
diff --git a/tests/scalapack/scalapack_10_b.with_hdf5=on.mpirun=1.output b/tests/scalapack/scalapack_10_b.with_hdf5=on.mpirun=1.output
new file mode 100644 (file)
index 0000000..ec9538e
--- /dev/null
@@ -0,0 +1,72 @@
+100x75 & 1 & 1x1
+100x75 & 1 & 10x10
+100x75 & 1 & 50x50
+100x75 & 1 & 100x75
+100x75 & 16 & 1x1
+100x75 & 16 & 10x10
+100x75 & 16 & 50x50
+100x75 & 16 & 100x75
+100x75 & 32 & 1x1
+100x75 & 32 & 10x10
+100x75 & 32 & 50x50
+100x75 & 32 & 100x75
+200x225 & 1 & 1x1
+200x225 & 1 & 10x10
+200x225 & 1 & 50x50
+200x225 & 1 & 100x75
+200x225 & 16 & 1x1
+200x225 & 16 & 10x10
+200x225 & 16 & 50x50
+200x225 & 16 & 100x75
+200x225 & 32 & 1x1
+200x225 & 32 & 10x10
+200x225 & 32 & 50x50
+200x225 & 32 & 100x75
+300x250 & 1 & 1x1
+300x250 & 1 & 10x10
+300x250 & 1 & 50x50
+300x250 & 1 & 100x75
+300x250 & 16 & 1x1
+300x250 & 16 & 10x10
+300x250 & 16 & 50x50
+300x250 & 16 & 100x75
+300x250 & 32 & 1x1
+300x250 & 32 & 10x10
+300x250 & 32 & 50x50
+300x250 & 32 & 100x75
+100x75 & 1 & 1x1
+100x75 & 1 & 10x10
+100x75 & 1 & 50x50
+100x75 & 1 & 100x75
+100x75 & 16 & 1x1
+100x75 & 16 & 10x10
+100x75 & 16 & 50x50
+100x75 & 16 & 100x75
+100x75 & 32 & 1x1
+100x75 & 32 & 10x10
+100x75 & 32 & 50x50
+100x75 & 32 & 100x75
+200x225 & 1 & 1x1
+200x225 & 1 & 10x10
+200x225 & 1 & 50x50
+200x225 & 1 & 100x75
+200x225 & 16 & 1x1
+200x225 & 16 & 10x10
+200x225 & 16 & 50x50
+200x225 & 16 & 100x75
+200x225 & 32 & 1x1
+200x225 & 32 & 10x10
+200x225 & 32 & 50x50
+200x225 & 32 & 100x75
+300x250 & 1 & 1x1
+300x250 & 1 & 10x10
+300x250 & 1 & 50x50
+300x250 & 1 & 100x75
+300x250 & 16 & 1x1
+300x250 & 16 & 10x10
+300x250 & 16 & 50x50
+300x250 & 16 & 100x75
+300x250 & 32 & 1x1
+300x250 & 32 & 10x10
+300x250 & 32 & 50x50
+300x250 & 32 & 100x75
diff --git a/tests/scalapack/scalapack_10_b.with_hdf5=on.mpirun=11.output b/tests/scalapack/scalapack_10_b.with_hdf5=on.mpirun=11.output
new file mode 100644 (file)
index 0000000..ec9538e
--- /dev/null
@@ -0,0 +1,72 @@
+100x75 & 1 & 1x1
+100x75 & 1 & 10x10
+100x75 & 1 & 50x50
+100x75 & 1 & 100x75
+100x75 & 16 & 1x1
+100x75 & 16 & 10x10
+100x75 & 16 & 50x50
+100x75 & 16 & 100x75
+100x75 & 32 & 1x1
+100x75 & 32 & 10x10
+100x75 & 32 & 50x50
+100x75 & 32 & 100x75
+200x225 & 1 & 1x1
+200x225 & 1 & 10x10
+200x225 & 1 & 50x50
+200x225 & 1 & 100x75
+200x225 & 16 & 1x1
+200x225 & 16 & 10x10
+200x225 & 16 & 50x50
+200x225 & 16 & 100x75
+200x225 & 32 & 1x1
+200x225 & 32 & 10x10
+200x225 & 32 & 50x50
+200x225 & 32 & 100x75
+300x250 & 1 & 1x1
+300x250 & 1 & 10x10
+300x250 & 1 & 50x50
+300x250 & 1 & 100x75
+300x250 & 16 & 1x1
+300x250 & 16 & 10x10
+300x250 & 16 & 50x50
+300x250 & 16 & 100x75
+300x250 & 32 & 1x1
+300x250 & 32 & 10x10
+300x250 & 32 & 50x50
+300x250 & 32 & 100x75
+100x75 & 1 & 1x1
+100x75 & 1 & 10x10
+100x75 & 1 & 50x50
+100x75 & 1 & 100x75
+100x75 & 16 & 1x1
+100x75 & 16 & 10x10
+100x75 & 16 & 50x50
+100x75 & 16 & 100x75
+100x75 & 32 & 1x1
+100x75 & 32 & 10x10
+100x75 & 32 & 50x50
+100x75 & 32 & 100x75
+200x225 & 1 & 1x1
+200x225 & 1 & 10x10
+200x225 & 1 & 50x50
+200x225 & 1 & 100x75
+200x225 & 16 & 1x1
+200x225 & 16 & 10x10
+200x225 & 16 & 50x50
+200x225 & 16 & 100x75
+200x225 & 32 & 1x1
+200x225 & 32 & 10x10
+200x225 & 32 & 50x50
+200x225 & 32 & 100x75
+300x250 & 1 & 1x1
+300x250 & 1 & 10x10
+300x250 & 1 & 50x50
+300x250 & 1 & 100x75
+300x250 & 16 & 1x1
+300x250 & 16 & 10x10
+300x250 & 16 & 50x50
+300x250 & 16 & 100x75
+300x250 & 32 & 1x1
+300x250 & 32 & 10x10
+300x250 & 32 & 50x50
+300x250 & 32 & 100x75
diff --git a/tests/scalapack/scalapack_10_b.with_hdf5=on.mpirun=4.output b/tests/scalapack/scalapack_10_b.with_hdf5=on.mpirun=4.output
new file mode 100644 (file)
index 0000000..ec9538e
--- /dev/null
@@ -0,0 +1,72 @@
+100x75 & 1 & 1x1
+100x75 & 1 & 10x10
+100x75 & 1 & 50x50
+100x75 & 1 & 100x75
+100x75 & 16 & 1x1
+100x75 & 16 & 10x10
+100x75 & 16 & 50x50
+100x75 & 16 & 100x75
+100x75 & 32 & 1x1
+100x75 & 32 & 10x10
+100x75 & 32 & 50x50
+100x75 & 32 & 100x75
+200x225 & 1 & 1x1
+200x225 & 1 & 10x10
+200x225 & 1 & 50x50
+200x225 & 1 & 100x75
+200x225 & 16 & 1x1
+200x225 & 16 & 10x10
+200x225 & 16 & 50x50
+200x225 & 16 & 100x75
+200x225 & 32 & 1x1
+200x225 & 32 & 10x10
+200x225 & 32 & 50x50
+200x225 & 32 & 100x75
+300x250 & 1 & 1x1
+300x250 & 1 & 10x10
+300x250 & 1 & 50x50
+300x250 & 1 & 100x75
+300x250 & 16 & 1x1
+300x250 & 16 & 10x10
+300x250 & 16 & 50x50
+300x250 & 16 & 100x75
+300x250 & 32 & 1x1
+300x250 & 32 & 10x10
+300x250 & 32 & 50x50
+300x250 & 32 & 100x75
+100x75 & 1 & 1x1
+100x75 & 1 & 10x10
+100x75 & 1 & 50x50
+100x75 & 1 & 100x75
+100x75 & 16 & 1x1
+100x75 & 16 & 10x10
+100x75 & 16 & 50x50
+100x75 & 16 & 100x75
+100x75 & 32 & 1x1
+100x75 & 32 & 10x10
+100x75 & 32 & 50x50
+100x75 & 32 & 100x75
+200x225 & 1 & 1x1
+200x225 & 1 & 10x10
+200x225 & 1 & 50x50
+200x225 & 1 & 100x75
+200x225 & 16 & 1x1
+200x225 & 16 & 10x10
+200x225 & 16 & 50x50
+200x225 & 16 & 100x75
+200x225 & 32 & 1x1
+200x225 & 32 & 10x10
+200x225 & 32 & 50x50
+200x225 & 32 & 100x75
+300x250 & 1 & 1x1
+300x250 & 1 & 10x10
+300x250 & 1 & 50x50
+300x250 & 1 & 100x75
+300x250 & 16 & 1x1
+300x250 & 16 & 10x10
+300x250 & 16 & 50x50
+300x250 & 16 & 100x75
+300x250 & 32 & 1x1
+300x250 & 32 & 10x10
+300x250 & 32 & 50x50
+300x250 & 32 & 100x75
diff --git a/tests/scalapack/scalapack_10_b.with_hdf5=on.mpirun=5.output b/tests/scalapack/scalapack_10_b.with_hdf5=on.mpirun=5.output
new file mode 100644 (file)
index 0000000..ec9538e
--- /dev/null
@@ -0,0 +1,72 @@
+100x75 & 1 & 1x1
+100x75 & 1 & 10x10
+100x75 & 1 & 50x50
+100x75 & 1 & 100x75
+100x75 & 16 & 1x1
+100x75 & 16 & 10x10
+100x75 & 16 & 50x50
+100x75 & 16 & 100x75
+100x75 & 32 & 1x1
+100x75 & 32 & 10x10
+100x75 & 32 & 50x50
+100x75 & 32 & 100x75
+200x225 & 1 & 1x1
+200x225 & 1 & 10x10
+200x225 & 1 & 50x50
+200x225 & 1 & 100x75
+200x225 & 16 & 1x1
+200x225 & 16 & 10x10
+200x225 & 16 & 50x50
+200x225 & 16 & 100x75
+200x225 & 32 & 1x1
+200x225 & 32 & 10x10
+200x225 & 32 & 50x50
+200x225 & 32 & 100x75
+300x250 & 1 & 1x1
+300x250 & 1 & 10x10
+300x250 & 1 & 50x50
+300x250 & 1 & 100x75
+300x250 & 16 & 1x1
+300x250 & 16 & 10x10
+300x250 & 16 & 50x50
+300x250 & 16 & 100x75
+300x250 & 32 & 1x1
+300x250 & 32 & 10x10
+300x250 & 32 & 50x50
+300x250 & 32 & 100x75
+100x75 & 1 & 1x1
+100x75 & 1 & 10x10
+100x75 & 1 & 50x50
+100x75 & 1 & 100x75
+100x75 & 16 & 1x1
+100x75 & 16 & 10x10
+100x75 & 16 & 50x50
+100x75 & 16 & 100x75
+100x75 & 32 & 1x1
+100x75 & 32 & 10x10
+100x75 & 32 & 50x50
+100x75 & 32 & 100x75
+200x225 & 1 & 1x1
+200x225 & 1 & 10x10
+200x225 & 1 & 50x50
+200x225 & 1 & 100x75
+200x225 & 16 & 1x1
+200x225 & 16 & 10x10
+200x225 & 16 & 50x50
+200x225 & 16 & 100x75
+200x225 & 32 & 1x1
+200x225 & 32 & 10x10
+200x225 & 32 & 50x50
+200x225 & 32 & 100x75
+300x250 & 1 & 1x1
+300x250 & 1 & 10x10
+300x250 & 1 & 50x50
+300x250 & 1 & 100x75
+300x250 & 16 & 1x1
+300x250 & 16 & 10x10
+300x250 & 16 & 50x50
+300x250 & 16 & 100x75
+300x250 & 32 & 1x1
+300x250 & 32 & 10x10
+300x250 & 32 & 50x50
+300x250 & 32 & 100x75

In the beginning the Universe was created. This has made a lot of people very angry and has been widely regarded as a bad move.

Douglas Adams


Typeset in Trocchi and Trocchi Bold Sans Serif.