def METHOD_NAME(refinement=1, ncycles=2, **kwargs):
    print_output(' --------\nRunning refinement {:}'.format(refinement))
    conservation_check = kwargs.pop('conservation_check', False)
    g_grav = float(physical_constants['g_grav'])
    depth = 100.0
    c_wave = numpy.sqrt(g_grav*depth)
    n_base = 20
    nx = n_base*refinement
    ny = 1
    lx = 60000.
    ly = lx/nx
    mesh2d = RectangleMesh(nx, ny, lx, ly)
    elev_amp = 10.0 if conservation_check else 0.01
    n_layers = 2*refinement
    # estimate of max advective velocity used to estimate time step
    u_mag = Constant(0.5)
    print_output('Triangle edge length {:} m'.format(lx/nx))
    print_output('Number of layers {:}'.format(n_layers))
    # bathymetry
    P1_2d = get_functionspace(mesh2d, 'CG', 1)
    bathymetry_2d = Function(P1_2d, name='Bathymetry')
    bathymetry_2d.assign(depth)
    # set time step, export interval and run duration
    n_steps = 1 if not conservation_check else 20
    T_cycle = lx/c_wave
    t_step = T_cycle/n_steps
    t_end = ncycles*T_cycle
    t_export = t_end if not conservation_check else t_step
    # create solver
    solver_obj = solver.FlowSolver(mesh2d, bathymetry_2d, n_layers)
    options = solver_obj.options
    options.element_family = 'dg-dg'
    options.timestepper_type = 'SSPRK22'
    options.use_nonlinear_equations = True
    options.solve_salinity = False
    options.solve_temperature = False
    options.use_implicit_vertical_diffusion = False
    options.use_baroclinic_formulation = False
    options.use_bottom_friction = False
    options.use_ale_moving_mesh = True
    options.use_limiter_for_tracers = True
    options.use_limiter_for_velocity = False
    options.simulation_export_time = t_export
    options.simulation_end_time = t_end
    options.horizontal_velocity_scale = u_mag
    options.check_volume_conservation_2d = conservation_check
    options.check_volume_conservation_3d = conservation_check
    options.check_salinity_conservation = conservation_check
    options.check_salinity_overshoot = conservation_check
    options.check_temperature_conservation = conservation_check
    options.check_temperature_overshoot = conservation_check
    options.fields_to_export = ['uv_2d', 'elev_2d', 'elev_3d', 'uv_3d',
                                'w_3d', 'w_mesh_3d', 'salt_3d', 'temp_3d',
                                'uv_dav_2d']
    options.fields_to_export_hdf5 = []
    options.update(kwargs)
    # need to call creator to create the function spaces
    solver_obj.create_equations()
    xy = SpatialCoordinate(solver_obj.mesh2d)
    elev_init = -elev_amp*cos(2*pi*xy[0]/lx)
    salt_init3d = temp_init3d = None
    if options.solve_salinity:
        salt_init3d = Function(solver_obj.function_spaces.H, name='initial salinity')
        salt_init3d.assign(4.5)
    if options.solve_temperature:
        temp_init3d = Function(solver_obj.function_spaces.H, name='initial temperature')
        x, y, z = SpatialCoordinate(solver_obj.mesh)
        temp_init3d.interpolate(5*sin(2*pi*x/lx) + 10.)
    solver_obj.assign_initial_conditions(elev=elev_init, salt=salt_init3d, temp=temp_init3d)
    solver_obj.iterate()
    area = lx*ly
    elev_err = errornorm(elev_init, solver_obj.fields.elev_2d)/numpy.sqrt(area)
    uv_err = errornorm(as_vector((0., 0.)), solver_obj.fields.uv_2d)/numpy.sqrt(area)
    print_output('L2 error elev={:}, uv={:}'.format(elev_err, uv_err))
    return elev_err, uv_err