anchor.diff:          end do
           call n_solvent_el%append(n_solvent)
   
  -        radial_hist%data = 0
  -        call compute_radial_histogram(radial_hist, colloids%pos(:,1), &
  -             solvent_cells%edges, solvent)
  -        call radial_hist_el%append(radial_hist%data)
  -
           v_com = sum(colloids%vel, dim=2)/colloids%Nmax
           call polar%update(com_pos, v_com, unit_r/norm2(unit_r), solvent, solvent_cells)
           call varia%tac()
  
positive.code: 
    type(h5md_element_t) :: q_el, omega_body_el, janus_pos_el, janus_vel_el, u_el
  
    type(histogram_t) :: radial_hist
    type(h5md_element_t) :: radial_hist_el
    double precision :: polar_r_max
    type(polar_fields_t) :: polar
    integer(HID_T) :: polar_id
  
negative.code:   end if
  
    call h5gcreate_f(hfile%id, 'fields', fields_group, error)
  
    call h5gcreate_f(hfile%id, 'connectivity', connectivity_group, error)
    call h5md_write_dataset(connectivity_group, 'janus_links', links-1)
  