canvases[s] = new CanvasCairo(scaling, universe->width(), universe->height());
   }
 
+  scalar_t gravity_fx = 0.0;
+  scalar_t gravity_fy = 1.0;
+
   scalar_t grab_start_x, grab_start_y;
 
   int failed;
       }
 
       if(s < nb_simulated_frames - 1) {
+
         // Run the simulation
+
         for(int i = 0; i < nb_iterations_per_steps; i++) {
           if(pulling) {
+            // Pulling the grabbed rectangle
             scalar_t xf = grabbed_polygon->absolute_x(grab_relative_x, grab_relative_y);
             scalar_t yf = grabbed_polygon->absolute_y(grab_relative_x, grab_relative_y);
             if (xf < 0 || xf >= world_width || yf < 0 || yf >= world_height) {
             grabbed_polygon->apply_force(dt, xf, yf, 0.0, -1.0);
           } else {
             // Gravity
+            universe->apply_gravity(dt, gravity_fx, gravity_fy);
           }
+
           universe->update(dt, 1.0 / scaling);
         }
       }
     total_nb_attempts++;
 
     if(total_nb_attempts >= max_total_nb_attempts) {
-      cerr << "There was " << max_total_nb_attempts << " attempts at generating the sequences." << endl;
+      cerr << "There was "
+           << max_total_nb_attempts
+           << " attempts at generating the sequences, aborting." << endl;
       abort();
     }
 
 
   }
 }
 
+void Universe::apply_gravity(scalar_t dt, scalar_t fx, scalar_t fy) {
+  for(int n = 0; n < _nb_polygons; n++)
+    if(_polygons[n])
+      _polygons[n]->apply_force(dt,
+                                _polygons[n]->_center_x, _polygons[n]->_center_y,
+                                fx, fy);
+}
+
 void Universe::apply_collision_forces(scalar_t dt) {
   const int nb_axis = 2;
   int nb_collision[_nb_polygons * _nb_polygons];
 
   // axis to speed up the computation
   void compute_pseudo_collisions(int nb_axis, int *nb_colliding_axis);
 
+  void apply_gravity(scalar_t dt, scalar_t fx, scalar_t fy);
   void apply_collision_forces(scalar_t dt);
   bool update(scalar_t dt, scalar_t padding);