anchor.diff:    double precision, dimension(3) :: center_P_initial, center_S_initial
   
     integer :: star_axis
  +  integer :: initial_motion_dir
   
     ! Density of ambient medium
     double precision :: ambient_density
  
positive.code:     use meth_params_module, only: rot_period
      use initial_model_module, only: init_1d
      use prob_params_module, only: center
  
      implicit none
  
      double precision :: dx
      type (eos_t) :: ambient_state
  
      call get_ambient(ambient_state)
  
  
negative.code: 
  
  
    ! Predict the locations of the primary and secondary stars at a given time assuming circular orbits.
  
    subroutine get_orbit_locations(time, loc_P, loc_S)
  
      use bl_error_module, only: bl_error
      use bl_constants_module, only: ZERO, TWO, FOUR, THIRD, M_PI
      use fundamental_constants_module, only: Gconst, M_solar
      use prob_params_module, only: center
      use meth_params_module, only: do_rotation, rot_period
  
      implicit none
  
      double precision, intent (in   ) :: time
      double precision, intent (inout) :: loc_P(3), loc_S(3)
  
      double precision :: theta
  
      ! Start off by assuming the stars are at the center of the grid.
      ! If we're only doing a single star simulation, then we can just stop here.
  
      loc_P = center
      loc_S = center
  
      if (mass_S < ZERO) return
  
      ! If rotation is enabled, then the stars should be at their initial locations.
      ! Otherwise, use the properties of a circular orbit to predict the locations.
  
      if (do_rotation .eq. 1) then
  
         theta = ZERO
  
       else
  
         theta = TWO * M_PI * time / rot_period
  
      endif
  
      loc_P(star_axis) = loc_P(star_axis) - a_P_initial * cos(theta)
      loc_S(star_axis) = loc_S(star_axis) + a_S_initial * cos(theta)
  
      loc_P(initial_motion_dir) = loc_P(initial_motion_dir) - a_P_initial * sin(theta)
      loc_S(initial_motion_dir) = loc_S(initial_motion_dir) + a_S_initial * sin(theta)
  
  
    end subroutine get_orbit_locations
  
  
  
    ! This routine checks to see if the primary mass is actually larger
    ! than the secondary mass, and switches them if not.
  
  