anchor.diff:    double precision :: r(3) !old direction vector
     double precision :: d
   
  -  d = sigma*2.d0 + 0.5d0!distance between the colloids in the dimer
  +  d = 2.d0*sigma + 0.5d0 !distance between the colloids in the dimer
     r_new = colloids% pos(:,1) - colloids% pos(:,2)
     s = colloids% pos_old(:,1) - colloids% pos_old(:,2)
     r = pos_rattle(:,1) - pos_rattle(:,2)
  
positive.code:   call solvent_cells%init(L, 1.d0)
    colloids% pos(:,1) = solvent_cells% edges/2.0
    colloids% pos(:,2) = solvent_cells% edges/2.0 
    colloids% pos(1,2) = colloids% pos(1,2) + 2.0*sigma + 0.5d0
    
    call solvent% random_placement(solvent_cells% edges, colloids, solvent_colloid_lj)
  
    call solvent% sort(solvent_cells)
  
    call neigh% init(colloids% Nmax, 300)
    call neigh% make_stencil(solvent_cells, 1.5d0)
  
    call neigh% update_list(colloids, solvent, 1.5d0, solvent_cells)
  
    tau = 0.1d0 !was 0.1d0 but we changed it to Kaprals' value 
    N_MD_steps = 100
  
negative.code:                  colloids% vel
  
       call solvent% sort(solvent_cells)
       call neigh% update_list(colloids, solvent, 1.5d0*sigma, solvent_cells)
  
       call simple_mpcd_step(solvent, solvent_cells, mt)
  
  