anchor.diff:  !
   !     end subroutine flux2
   
  -    subroutine advect3d(q,u,v,w,rho,dz,nx,nz,ny,options)
  +    subroutine advect3d(q,u,v,w,rho,dz,dx,nx,nz,ny,options)
           implicit none
           real, dimension(1:nx,  1:nz,1:ny),  intent(inout)   :: q
           real, dimension(1:nx,  1:nz,1:ny),  intent(in)      :: w
           real, dimension(1:nx-1,1:nz,1:ny),  intent(in)      :: u
           real, dimension(1:nx,  1:nz,1:ny-1),intent(in)      :: v
  -        real, dimension(:, :, :),  intent(in)      :: rho
  -        real, dimension(:, :, :),  intent(in)      :: dz
  +        real, dimension(:, :, :),           intent(in)      :: rho
  +        real, dimension(:, :, :),           intent(in)      :: dz
  +        real,                               intent(in)      :: dx
           integer,                            intent(in)      :: ny, nz, nx
           type(options_t), intent(in)::options
   
  
positive.code:     !
      ! end subroutine advect_cu_winds
  
      subroutine setup_module_winds(u,v,w, dx, options, dt)
          implicit none
          real,               intent(in)  :: u(:,:,:)
          real,               intent(in)  :: v(:,:,:)
          real,               intent(in)  :: w(:,:,:)
          real,               intent(in)  :: dx
          type(options_t),    intent(in)  :: options
          real,               intent(in)  :: dt
  
negative.code:         ny = domain%grid%ny
  
          ! calculate U,V,W normalized for dt/dx (dx**2 for density advection so we can skip a /dx in the actual advection code)
          call setup_module_winds(domain%u%data_3d, domain%v%data_3d, domain%w%data_3d, domain%dz_interface%data_3d, domain%dx, options, dt)
  
          ! lastqv_m=domain%qv
  
          if (options%vars_to_advect(kVARS%water_vapor)>0)                  call advect3d(domain%water_vapor%data_3d,             U_m,V_m,W_m, domain%density%data_3d, domain%dz_interface%data_3d, domain%dx, nx,nz,ny, options)
          if (options%vars_to_advect(kVARS%cloud_water)>0)                  call advect3d(domain%cloud_water_mass%data_3d,        U_m,V_m,W_m, domain%density%data_3d, domain%dz_interface%data_3d, domain%dx, nx,nz,ny, options)
          if (options%vars_to_advect(kVARS%rain_in_air)>0)                  call advect3d(domain%rain_mass%data_3d,               U_m,V_m,W_m, domain%density%data_3d, domain%dz_interface%data_3d, domain%dx, nx,nz,ny, options)
          if (options%vars_to_advect(kVARS%snow_in_air)>0)                  call advect3d(domain%snow_mass%data_3d,               U_m,V_m,W_m, domain%density%data_3d, domain%dz_interface%data_3d, domain%dx, nx,nz,ny, options)
          if (options%vars_to_advect(kVARS%potential_temperature)>0)        call advect3d(domain%potential_temperature%data_3d,   U_m,V_m,W_m, domain%density%data_3d, domain%dz_interface%data_3d, domain%dx, nx,nz,ny, options)
          if (options%vars_to_advect(kVARS%cloud_ice)>0)                    call advect3d(domain%cloud_ice_mass%data_3d,          U_m,V_m,W_m, domain%density%data_3d, domain%dz_interface%data_3d, domain%dx, nx,nz,ny, options)
          if (options%vars_to_advect(kVARS%graupel_in_air)>0)               call advect3d(domain%graupel_mass%data_3d,            U_m,V_m,W_m, domain%density%data_3d, domain%dz_interface%data_3d, domain%dx, nx,nz,ny, options)
          if (options%vars_to_advect(kVARS%ice_number_concentration)>0)     call advect3d(domain%cloud_ice_number%data_3d,        U_m,V_m,W_m, domain%density%data_3d, domain%dz_interface%data_3d, domain%dx, nx,nz,ny, options)
          if (options%vars_to_advect(kVARS%rain_number_concentration)>0)    call advect3d(domain%rain_number%data_3d,             U_m,V_m,W_m, domain%density%data_3d, domain%dz_interface%data_3d, domain%dx, nx,nz,ny, options)
          if (options%vars_to_advect(kVARS%snow_number_concentration)>0)    call advect3d(domain%snow_number%data_3d,             U_m,V_m,W_m, domain%density%data_3d, domain%dz_interface%data_3d, domain%dx, nx,nz,ny, options)
          if (options%vars_to_advect(kVARS%graupel_number_concentration)>0) call advect3d(domain%graupel_number%data_3d,          U_m,V_m,W_m, domain%density%data_3d, domain%dz_interface%data_3d, domain%dx, nx,nz,ny, options)
  
          ! if (options%physics%convection > 0) then
          !     call advect_cu_winds(domain, options, dt)
  