anchor.diff:  
     tau = PTread_d(config, 'tau', loc=params_group)
     N_MD_steps = PTread_i(config, 'N_MD', loc=params_group)
  +  colloid_sampling = PTread_i(config, 'colloid_sampling', loc=params_group)
  +  if (modulo(N_MD_steps, colloid_sampling) /= 0) then
  +     error stop 'colloid_sampling must divide N_MD with no remainder'
  +  end if
     dt = tau / N_MD_steps
     N_loop = PTread_i(config, 'N_loop', loc=params_group)
     steps_fixed = PTread_i(config, 'steps_fixed', loc=params_group)
  
positive.code:   dimer_io%id_info%store = .false.
    dimer_io%position_info%store = .true.
    dimer_io%position_info%mode = ior(H5MD_LINEAR,H5MD_STORE_TIME)
    dimer_io%position_info%step = N_MD_steps
    dimer_io%position_info%time = N_MD_steps*dt
    dimer_io%image_info%store = .true.
    dimer_io%image_info%mode = ior(H5MD_LINEAR,H5MD_STORE_TIME)
    dimer_io%image_info%step = N_MD_steps
    dimer_io%image_info%time = N_MD_steps*dt
    dimer_io%velocity_info%store = .true.
    dimer_io%velocity_info%mode = ior(H5MD_LINEAR,H5MD_STORE_TIME)
    dimer_io%velocity_info%step = N_MD_steps
    dimer_io%velocity_info%time = N_MD_steps*dt
    dimer_io%species_info%store = .true.
    dimer_io%species_info%mode = H5MD_FIXED
    call dimer_io%init(hfile, 'dimer', colloids)
  
negative.code:            end if
          end if
  
          if (modulo(j, colloid_sampling) == 0) then
             call dimer_io%position%append(colloids%pos)
             call dimer_io%velocity%append(colloids%vel)
             call dimer_io%image%append(colloids%image)
          end if
  
       end do md_loop
  
       call random_number(solvent_cells% origin)
  