anchor.diff:  
     call main%init('main')
     call varia%init('varia')
  +  call so_timer%init('so_timer')
     call time_flag%init('flag')
     call time_refuel%init('refuel')
     call time_change%init('change')
     call bulk_reac_timer%init('bulk_reac')
     call q_timer%init('quaternion_vv')
   
  -  call timer_list%init(23)
  +  call timer_list%init(24)
     call timer_list%append(varia)
  +  call timer_list%append(so_timer)
     call timer_list%append(time_flag)
     call timer_list%append(time_refuel)
     call timer_list%append(time_change)
  
positive.code:            call apply_pbc(colloids, solvent_cells% edges)
             call solvent% sort(solvent_cells)
             call neigh% update_list(colloids, solvent, max_cut + skin, solvent_cells, solvent_colloid_lj)
             call varia%tic()
             !$omp parallel do
             do k = 1, solvent%Nmax
                solvent% pos_old(:,k) = solvent% pos(:,k)
             end do
             colloids% pos_old = colloids% pos
             n_extra_sorting = n_extra_sorting + 1
             call varia%tac()
  
negative.code:         call varia%tic()
          call switch(solvent% force, solvent% force_old)
          call switch(colloids% force, colloids% force_old)
          call varia%tac()
  
          call so_timer%tic()
          !$omp parallel do
          do k = 1, solvent%Nmax
             solvent% force(:,k) = 0
          end do
          call so_timer%tac()
          colloids% force = 0
  
          e1 = compute_force(colloids, solvent, neigh, solvent_cells% edges, solvent_colloid_lj)
          if (do_lennard_jones) e2 = compute_force_n2(colloids, solvent_cells% edges, colloid_lj)
          if (do_elastic) e3 = elastic_network(colloids, links, links_d, elastic_k, solvent_cells%edges)
  