A multi-mode nonlinear constitutive model for mucus is constructed directly from micro- and macro-rheology experimental data on cell culture mucus, and a numerical algorithm is developed for the culture geometry and idealized cilia driving conditions. This study investigates the roles that mucus rheology, wall effects, and HBE culture geometry play in the development of flow profiles and the shape of the air-mucus interface. Simulations show that viscoelasticity captures normal stress generation in shear leading to a peak in the air-mucus interface at the middle of the culture and a depression at the walls. Linear and nonlinear viscoelastic regimes can be observed in cultures by varying the hurricane radius and mean rotational velocity. The advection-diffusion of a drug concentration dropped at the surface of the mucus flow is simulated as a function of Peclet number.