]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
Shuffle around a few functions that can't be shared between sequential and parallel...
authorwolf <wolf@0785d39b-7218-0410-832d-ea1e28bc413d>
Wed, 24 Mar 2004 19:38:53 +0000 (19:38 +0000)
committerwolf <wolf@0785d39b-7218-0410-832d-ea1e28bc413d>
Wed, 24 Mar 2004 19:38:53 +0000 (19:38 +0000)
git-svn-id: https://svn.dealii.org/trunk@8861 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/lac/include/lac/petsc_parallel_vector.h
deal.II/lac/include/lac/petsc_vector.h
deal.II/lac/include/lac/petsc_vector_base.h
deal.II/lac/source/petsc_parallel_vector.cc
deal.II/lac/source/petsc_vector.cc
deal.II/lac/source/petsc_vector_base.cc

index e7c12e7453bedd28ee5aaf5a5c0569877cded550..8f13f6c3b9f3ff9c9edaf58cb583f375d152fa29 100644 (file)
@@ -61,13 +61,14 @@ namespace PETScWrappers
                                           * @p{n} and initialize all
                                           * elements with zero.
                                           *
-                                          * The constructor is made explicit to
-                                          * avoid accidents like this:
-                                          * @p{v=0;}. Presumably, the user wants
-                                          * to set every element of the vector to
-                                          * zero, but instead, what happens is
-                                          * this call: @p{v=Vector<number>(0);},
-                                          * i.e. the vector is replaced by one of
+                                          * The constructor is made explicit
+                                          * to avoid accidents like this:
+                                          * @p{v=0;}. Presumably, the user
+                                          * wants to set every element of the
+                                          * vector to zero, but instead, what
+                                          * happens is this call:
+                                          * @p{v=Vector<number>(0);}, i.e. the
+                                          * vector is replaced by one of
                                           * length zero.
                                           */
         explicit Vector (const unsigned int n,
@@ -94,6 +95,14 @@ namespace PETScWrappers
                          const unsigned int local_size,
                          const MPI_Comm    &communicator);
 
+                                         /**
+                                          * Copy the given vector. Resize the
+                                          * present vector if necessary. Also
+                                          * take over the MPI communicator of
+                                          * @arg v.
+                                          */
+        Vector & operator = (const Vector &v);
+
                                          /**
                                           * Set all components of the vector to
                                           * the given number @p{s}. Simply pass
@@ -110,20 +119,76 @@ namespace PETScWrappers
                                           * (as opposed to those of the PETSc
                                           * vector wrapper class) into this
                                           * object.
+                                          *
+                                          * Contrary to the case of sequential
+                                          * vectors, this operators requires
+                                          * that the present vector already
+                                          * has the correct size, since we
+                                          * need to have a partition and a
+                                          * communicator present which we
+                                          * otherwise can't get from the
+                                          * source vector.
                                           */
         template <typename number>
         Vector & operator = (const ::Vector<number> &v);
-      
+
+                                         /**
+                                          * Change the dimension of the vector
+                                          * to @arg N. It is unspecified how
+                                          * resizing the vector affects the
+                                          * memory allocation of this object;
+                                          * i.e., it is not guaranteed that
+                                          * resizing it to a smaller size
+                                          * actually also reduces memory
+                                          * consumption, or if for efficiency
+                                          * the same amount of memory is used
+                                          *
+                                          * @arg local_size denotes how many
+                                          * of the @arg N values shall be
+                                          * stored locally on the present
+                                          * process.
+                                          * for less data.
+                                          *
+                                          * @arg communicator denotes the MPI
+                                          * communicator henceforth to be used
+                                          * for this vector.
+                                          * 
+                                          * If @arg fast is false, the vector
+                                          * is filled by zeros. Otherwise, the
+                                          * elements are left an unspecified
+                                          * state.
+                                          */ 
+        void reinit (const unsigned int N,
+                     const unsigned int local_size,
+                     const MPI_Comm    &communicator,
+                     const bool         fast = false);
+    
+                                         /**
+                                          * Change the dimension to that of
+                                          * the vector @arg v, and also take
+                                          * over the partitioning into local
+                                          * sizes as well as the MPI
+                                          * communicator. The same applies as
+                                          * for the other @p{reinit} function.
+                                          *
+                                          * The elements of @arg v are not
+                                          * copied, i.e. this function is the
+                                          * same as calling
+                                          * <tt>reinit(v.size(),
+                                          * v.local_size(), fast)</tt>.
+                                          */
+        void reinit (const Vector &v,
+                     const bool    fast = false);
+
       protected:
                                          /**
-                                          * Create a vector of length @p{n}. For
-                                          * this class, we create a parallel
-                                          * vector. @arg n denotes the total
-                                          * size of the vector to be
-                                          * created. @arg local_size denotes how
-                                          * many of these elements shall be
-                                          * stored locally. The last argument is
-                                          * ignored for sequential vectors.
+                                          * Create a vector of length
+                                          * @p{n}. For this class, we create a
+                                          * parallel vector. @arg n denotes
+                                          * the total size of the vector to be
+                                          * created. @arg local_size denotes
+                                          * how many of these elements shall
+                                          * be stored locally.
                                           */
         virtual void create_vector (const unsigned int n,
                                     const unsigned int local_size);
@@ -150,7 +215,7 @@ namespace PETScWrappers
     {
       Vector::create_vector (v.size(), local_size);
 
-      VectorBase::operator = (v);
+      *this = v;
     }
 
   
@@ -165,16 +230,81 @@ namespace PETScWrappers
     }
   
 
+
+    inline
+    Vector &
+    Vector::operator = (const Vector &v)
+    {
+                                       // if the vectors have different sizes,
+                                       // then first resize the present one
+      if (size() != v.size())
+        reinit (v.size(), v.local_size(), v.communicator, true);
+    
+      const int ierr = VecCopy (v.vector, vector);
+      AssertThrow (ierr == 0, ExcPETScError(ierr));
+    
+      return *this;
+    }
+
+
+
     template <typename number>
     inline
     Vector &
-    Vector::operator = (const ::Vector<number> &v)
+    Vector::operator = (const ::Vector<number> &v) 
     {
-      VectorBase::operator = (v);
+      Assert (size() == v.size(),
+              ExcNonMatchingSizes (size(), v.size()));
+
+                                       // the following isn't necessarily fast,
+                                       // but this is due to the fact that PETSc
+                                       // doesn't offer an inlined access
+                                       // operator.
+                                       //
+                                       // if someone wants to contribute some
+                                       // code: to make this code faster, one
+                                       // could either first convert all values
+                                       // to PetscScalar, and then set them all
+                                       // at once using VecSetValues. This has
+                                       // the drawback that it could take quite
+                                       // some memory, if the vector is large,
+                                       // and it would in addition allocate
+                                       // memory on the heap, which is
+                                       // expensive. an alternative would be to
+                                       // split the vector into chunks of, say,
+                                       // 128 elements, convert a chunk at a
+                                       // time and set it in the output vector
+                                       // using VecSetValues. since 128 elements
+                                       // is small enough, this could easily be
+                                       // allocated on the stack (as a local
+                                       // variable) which would make the whole
+                                       // thing much more efficient.
+                                       //
+                                       // a second way to make things faster is
+                                       // for the special case that
+                                       // number==PetscScalar. we could then
+                                       // declare a specialization of this
+                                       // template, and omit the conversion. the
+                                       // problem with this is that the best we
+                                       // can do is to use VecSetValues, but
+                                       // this isn't very efficient either: it
+                                       // wants to see an array of indices,
+                                       // which in this case a) again takes up a
+                                       // whole lot of memory on the heap, and
+                                       // b) is totally dumb since its content
+                                       // would simply be the sequence
+                                       // 0,1,2,3,...,n. the best of all worlds
+                                       // would probably be a function in Petsc
+                                       // that would take a pointer to an array
+                                       // of PetscScalar values and simply copy
+                                       // n elements verbatim into the vector...
+      for (unsigned int i=0; i<v.size(); ++i)
+        (*this)(i) = v(i);
+
+      compress ();
 
       return *this;
-    }
-    
+    }    
   }
   
 }
index f23a6bae4895848dc0ff95cd936d3709d5629dd8..f9786c918c046bad4d38284d04a322b524f01b98 100644 (file)
@@ -78,6 +78,12 @@ namespace PETScWrappers
                                         */
       explicit Vector (const VectorBase &v);
 
+                                       /**
+                                        * Copy the given vector. Resize the
+                                        * present vector if necessary.
+                                        */
+      Vector & operator = (const Vector &v);
+
                                        /**
                                         * Set all components of the vector to
                                         * the given number @p{s}. Simply pass
@@ -97,20 +103,49 @@ namespace PETScWrappers
                                         */
       template <typename number>
       Vector & operator = (const ::Vector<number> &v);
-      
+
+                                       /**
+                                        * Change the dimension of the vector
+                                        * to @arg N. It is unspecified how
+                                        * resizing the vector affects the
+                                        * memory allocation of this object;
+                                        * i.e., it is not guaranteed that
+                                        * resizing it to a smaller size
+                                        * actually also reduces memory
+                                        * consumption, or if for efficiency
+                                        * the same amount of memory is used
+                                        * for less data.
+                                        *
+                                        * If @arg fast is false, the vector is
+                                        * filled by zeros. Otherwise, the
+                                        * elements are left an unspecified
+                                        * state.
+                                        */ 
+      void reinit (const unsigned int N,
+                   const bool         fast = false);
+    
+                                       /**
+                                        * Change the dimension to that of the
+                                        * vector @arg v. The same applies as
+                                        * for the other reinit() function.
+                                        *
+                                        * The elements of @arg v are not
+                                        * copied, i.e.  this function is the
+                                        * same as calling <tt>reinit (v.size(),
+                                        * fast)</tt>.
+                                        */
+      void reinit (const Vector &v,
+                   const bool    fast = false);
+
     protected:
                                        /**
                                         * Create a vector of length @p{n}. For
-                                        * this class, we create a parallel
+                                        * this class, we create a sequential
                                         * vector. @arg n denotes the total
                                         * size of the vector to be
-                                        * created. @arg local_size denotes how
-                                        * many of these elements shall be
-                                        * stored locally. The last argument is
-                                        * ignored for sequential vectors.
+                                        * created.
                                         */
-      virtual void create_vector (const unsigned int n,
-                                  const unsigned int local_size = 0);
+      void create_vector (const unsigned int n);
   };
 
 
@@ -123,7 +158,7 @@ namespace PETScWrappers
   {
     Vector::create_vector (v.size());
 
-    VectorBase::operator = (v);
+    *this = v;
   }
 
   
@@ -138,16 +173,78 @@ namespace PETScWrappers
   }
   
 
+  inline
+  Vector &
+  Vector::operator = (const Vector &v)
+  {
+                                     // if the vectors have different sizes,
+                                     // then first resize the present one
+    if (size() != v.size())
+      reinit (v.size(), true);
+    
+    const int ierr = VecCopy (v.vector, vector);
+    AssertThrow (ierr == 0, ExcPETScError(ierr));
+    
+    return *this;
+  }
+
+
+
   template <typename number>
   inline
   Vector &
-  Vector::operator = (const ::Vector<number> &v)
+  Vector::operator = (const ::Vector<number> &v) 
   {
-    VectorBase::operator = (v);
+    reinit (v.size(), true);
+                                     // the following isn't necessarily fast,
+                                     // but this is due to the fact that PETSc
+                                     // doesn't offer an inlined access
+                                     // operator.
+                                     //
+                                     // if someone wants to contribute some
+                                     // code: to make this code faster, one
+                                     // could either first convert all values
+                                     // to PetscScalar, and then set them all
+                                     // at once using VecSetValues. This has
+                                     // the drawback that it could take quite
+                                     // some memory, if the vector is large,
+                                     // and it would in addition allocate
+                                     // memory on the heap, which is
+                                     // expensive. an alternative would be to
+                                     // split the vector into chunks of, say,
+                                     // 128 elements, convert a chunk at a
+                                     // time and set it in the output vector
+                                     // using VecSetValues. since 128 elements
+                                     // is small enough, this could easily be
+                                     // allocated on the stack (as a local
+                                     // variable) which would make the whole
+                                     // thing much more efficient.
+                                     //
+                                     // a second way to make things faster is
+                                     // for the special case that
+                                     // number==PetscScalar. we could then
+                                     // declare a specialization of this
+                                     // template, and omit the conversion. the
+                                     // problem with this is that the best we
+                                     // can do is to use VecSetValues, but
+                                     // this isn't very efficient either: it
+                                     // wants to see an array of indices,
+                                     // which in this case a) again takes up a
+                                     // whole lot of memory on the heap, and
+                                     // b) is totally dumb since its content
+                                     // would simply be the sequence
+                                     // 0,1,2,3,...,n. the best of all worlds
+                                     // would probably be a function in Petsc
+                                     // that would take a pointer to an array
+                                     // of PetscScalar values and simply copy
+                                     // n elements verbatim into the vector...
+    for (unsigned int i=0; i<v.size(); ++i)
+      (*this)(i) = v(i);
+
+    compress ();
 
     return *this;
   }
-  
 }
 
 
index 6dae8e473518b7154fa11a95e0270700721e85e5..796d1d71714ccad9f77108aa8103b85b7c20b79b 100644 (file)
@@ -208,47 +208,6 @@ namespace PETScWrappers
                                         */
       void compress ();
 
-                                       /**
-                                        * Change the dimension of the vector
-                                        * to @arg N. It is unspecified how
-                                        * resizing the vector affects the
-                                        * memory allocation of this object;
-                                        * i.e., it is not guaranteed that
-                                        * resizing it to a smaller size
-                                        * actually also reduces memory
-                                        * consumption, or if for efficiency
-                                        * the same amount of memory is used
-                                        * for less data.
-                                        *
-                                        * On @arg fast is false, the vector is
-                                        * filled by zeros. Otherwise, the
-                                        * elements are left an unspecified
-                                        * state.
-                                        *
-                                        * For parallel vectors, @arg
-                                        * local_size denotes how many of the
-                                        * @arg N values shall be stored
-                                        * locally on the present process. This
-                                        * argument is ignored for sequantial
-                                        * vectors.
-                                        */ 
-      void reinit (const unsigned int N,
-                   const bool         fast = false,
-                   const unsigned int local_size = 0);
-    
-                                       /**
-                                        * Change the dimension to that of the
-                                        * vector @p{V}. The same applies as
-                                        * for the other @p{reinit} function.
-                                        *
-                                        * The elements of @p{V} are not
-                                        * copied, i.e.  this function is the
-                                        * same as calling @p{reinit (V.size(),
-                                        * fast)}.
-                                        */
-      void reinit (const VectorBase &V,
-                   const bool        fast = false);
-
                                        /**
                                         * Set all entries to zero. Equivalent
                                         * to @p{v = 0}, but more obvious and
@@ -263,22 +222,7 @@ namespace PETScWrappers
                                         * Set all components of the vector to
                                         * the given number @p{s}.
                                         */
-      VectorBase & operator = (const PetscScalar s);
-    
-                                       /**
-                                        * Copy the given vector. Resize the
-                                        * present vector if necessary.
-                                        */
-      VectorBase & operator = (const VectorBase &v);
-
-                                       /**
-                                        * Copy the values of a deal.II vector
-                                        * (as opposed to those of the PETSc
-                                        * vector wrapper class) into this
-                                        * object.
-                                        */
-      template <typename number>
-      VectorBase & operator = (const ::Vector<number> &v);
+      VectorBase & operator = (const PetscScalar s);    
       
                                        /**
                                         * Test for equality. This function
@@ -588,7 +532,14 @@ namespace PETScWrappers
                       int,
                       << "An error with error number " << arg1
                       << " occured while calling a PETSc function");
-
+                                       /**
+                                        * Exception
+                                        */
+      DeclException2 (ExcNonMatchingSizes,
+                      int, int,
+                      << "The sizes " << arg1 << " and " << arg2
+                      << " are supposed to be equal, but are not.");
+      
     protected:
                                        /**
                                         * A generic vector object in
@@ -635,19 +586,6 @@ namespace PETScWrappers
                                         */
       mutable LastAction::Values last_action;
 
-                                       /**
-                                        * Create a vector of length @p{n}. For
-                                        * this class, we create a parallel
-                                        * vector. @arg n denotes the total
-                                        * size of the vector to be
-                                        * created. @arg local_size denotes how
-                                        * many of these elements shall be
-                                        * stored locally. The last argument is
-                                        * ignored for sequential vectors.
-                                        */
-      virtual void create_vector (const unsigned int n,
-                                  const unsigned int local_size = 0) = 0;
-
                                        /**
                                         * Make the reference class a friend.
                                         */
@@ -876,60 +814,6 @@ namespace PETScWrappers
   
 
 
-  template <typename number>
-  VectorBase &
-  VectorBase::operator = (const ::Vector<number> &v) 
-  {
-    reinit (v.size());
-                                     // the following isn't necessarily fast,
-                                     // but this is due to the fact that PETSc
-                                     // doesn't offer an inlined access
-                                     // operator.
-                                     //
-                                     // if someone wants to contribute some
-                                     // code: to make this code faster, one
-                                     // could either first convert all values
-                                     // to PetscScalar, and then set them all
-                                     // at once using VecSetValues. This has
-                                     // the drawback that it could take quite
-                                     // some memory, if the vector is large,
-                                     // and it would in addition allocate
-                                     // memory on the heap, which is
-                                     // expensive. an alternative would be to
-                                     // split the vector into chunks of, say,
-                                     // 128 elements, convert a chunk at a
-                                     // time and set it in the output vector
-                                     // using VecSetValues. since 128 elements
-                                     // is small enough, this could easily be
-                                     // allocated on the stack (as a local
-                                     // variable) which would make the whole
-                                     // thing much more efficient.
-                                     //
-                                     // a second way to make things faster is
-                                     // for the special case that
-                                     // number==PetscScalar. we could then
-                                     // declare a specialization of this
-                                     // template, and omit the conversion. the
-                                     // problem with this is that the best we
-                                     // can do is to use VecSetValues, but
-                                     // this isn't very efficient either: it
-                                     // wants to see an array of indices,
-                                     // which in this case a) again takes up a
-                                     // whole lot of memory on the heap, and
-                                     // b) is totally dumb since its content
-                                     // would simply be the sequence
-                                     // 0,1,2,3,...,n. the best of all worlds
-                                     // would probably be a function in Petsc
-                                     // that would take a pointer to an array
-                                     // of PetscScalar values and simply copy
-                                     // n elements verbatim into the vector...
-    for (unsigned int i=0; i<v.size(); ++i)
-      (*this)(i) = v(i);
-
-    compress ();
-
-    return *this;
-  }
 }
 
 #endif // DEAL_II_USE_PETSC
index 93371a56891cf8c2437898c2c47b3cd0558cc975..9d6b77abba3d8a4e4704022ab182bea6b5a5f95b 100644 (file)
@@ -61,6 +61,53 @@ namespace PETScWrappers
     }
 
   
+
+    void
+    Vector::reinit (const unsigned int n,
+                    const unsigned int local_sz,
+                    const MPI_Comm    &comm,
+                    const bool         fast)
+    {
+      communicator = comm;
+      
+                                       // only do something if the sizes
+                                       // mismatch
+      if ((size() != n) || (local_size() != local_sz))
+        {
+                                           // FIXME: I'd like to use this here,
+                                           // but somehow it leads to odd errors
+                                           // somewhere down the line in some of
+                                           // the tests:
+//         const int ierr = VecSetSizes (vector, n, n);
+//         AssertThrow (ierr == 0, ExcPETScError(ierr));
+
+                                           // so let's go the slow way:
+          int ierr;
+          ierr = VecDestroy (vector);
+          AssertThrow (ierr == 0, ExcPETScError(ierr));
+
+          create_vector (n, local_sz);
+        }
+
+                                       // finally clear the new vector if so
+                                       // desired
+      if (fast == false)
+        *this = 0;
+    }
+
+
+
+    void
+    Vector::reinit (const Vector &v,
+                    const bool    fast)
+    {
+      communicator = v.communicator;
+      
+      reinit (v.size(), v.local_size(), fast);
+    }
+  
+
+
     void
     Vector::create_vector (const unsigned int  n,
                            const unsigned int  local_size)
index 43484710b110c9d32e71aacd628cb11b891e18ef..29633938ddf7386d910d2ab37be290c31e64db21 100644 (file)
@@ -43,13 +43,51 @@ namespace PETScWrappers
     VectorBase::operator = (v);
   }
 
+
+
+  void
+  Vector::reinit (const unsigned int n,
+                  const bool         fast)
+  {
+                                     // only do something if the sizes
+                                     // mismatch
+    if (size() != n)
+      {
+                                         // FIXME: I'd like to use this here,
+                                         // but somehow it leads to odd errors
+                                         // somewhere down the line in some of
+                                         // the tests:
+//         const int ierr = VecSetSizes (vector, n, n);
+//         AssertThrow (ierr == 0, ExcPETScError(ierr));
+
+                                         // so let's go the slow way:
+        int ierr;
+        ierr = VecDestroy (vector);
+        AssertThrow (ierr == 0, ExcPETScError(ierr));
+
+        create_vector (n);
+      }
+
+                                     // finally clear the new vector if so
+                                     // desired
+    if (fast == false)
+      *this = 0;
+  }
+
+
+
+  void
+  Vector::reinit (const Vector &v,
+                  const bool    fast)
+  {
+    reinit (v.size(), fast);
+  }
+  
+
   
   void
-  Vector::create_vector (const unsigned int n,
-                         const unsigned int local_size)
+  Vector::create_vector (const unsigned int n)
   {
-    Assert (local_size < n, ExcIndexRange (local_size, 0, n));
-    
     const int ierr
       = VecCreateSeq (PETSC_COMM_SELF, n, &vector);
     AssertThrow (ierr == 0, ExcPETScError(ierr));
index b074a01fdbb8072ba6521fda42efd9d98e908226..f251402640ce5514eecb1d6f92ff625c5226bcb6 100644 (file)
@@ -48,47 +48,6 @@ namespace PETScWrappers
 
 
 
-  void
-  VectorBase::reinit (const unsigned int n,
-                      const bool         fast,
-                      const unsigned int local_sz)
-  {
-                                     // only do something if the sizes
-                                     // mismatch
-    if ((size() != n) || (local_size() != local_sz))
-      {
-                                         // FIXME: I'd like to use this here,
-                                         // but somehow it leads to odd errors
-                                         // somewhere down the line in some of
-                                         // the tests:
-//         const int ierr = VecSetSizes (vector, n, n);
-//         AssertThrow (ierr == 0, ExcPETScError(ierr));
-
-                                         // so let's go the slow way:
-        int ierr;
-        ierr = VecDestroy (vector);
-        AssertThrow (ierr == 0, ExcPETScError(ierr));
-
-        create_vector (n, local_sz);
-      }
-
-                                     // finally clear the new vector if so
-                                     // desired
-    if (fast == false)
-      *this = 0;
-  }
-
-
-
-  void
-  VectorBase::reinit (const VectorBase &v,
-                      const bool        fast)
-  {
-    reinit (v.size(), fast);
-  }
-  
-
-  
   void
   VectorBase::clear ()
   {
@@ -110,22 +69,6 @@ namespace PETScWrappers
 
 
 
-  VectorBase &
-  VectorBase::operator = (const VectorBase &v)
-  {
-                                     // if the vectors have different sizes,
-                                     // then first resize the present one
-    if (size() != v.size())
-      reinit (v.size(), true);
-    
-    const int ierr = VecCopy (v.vector, vector);
-    AssertThrow (ierr == 0, ExcPETScError(ierr));
-
-    return *this;
-  }
-  
-
-
   bool
   VectorBase::operator == (const VectorBase &v) const
   {
@@ -595,10 +538,15 @@ namespace PETScWrappers
   VectorBase::equ (const PetscScalar a,
                    const VectorBase     &v)
   {
+    Assert (size() == v.size(),
+            ExcNonMatchingSizes (size(), v.size()));
+
                                      // there is no simple operation for this
                                      // in PETSc. there are multiple ways to
                                      // emulate it, we choose this one:
-    *this = v;
+    const int ierr = VecCopy (v.vector, vector);
+    AssertThrow (ierr == 0, ExcPETScError(ierr));
+
     *this *= a;
   }
   
@@ -610,10 +558,15 @@ namespace PETScWrappers
                    const PetscScalar b,
                    const VectorBase     &w)
   {
+    Assert (size() == v.size(),
+            ExcNonMatchingSizes (size(), v.size()));
+
                                      // there is no simple operation for this
                                      // in PETSc. there are multiple ways to
                                      // emulate it, we choose this one:
-    *this = v;
+    const int ierr = VecCopy (v.vector, vector);
+    AssertThrow (ierr == 0, ExcPETScError(ierr));
+
     sadd (a, b, w);
   }
 

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.