This paper makes use of sloshing reduced-order models to investigate the effects of sloshing dynamics on aeroelastic stability and response of flying wing structure. More specifically, a linear frequency-domain operator derived by an equivalent mechanical model is used to model lateral (linear) sloshing dynamics whereas data-driven neural-networks are used to model the vertical (nonlinear) sloshing dynamics. These models are integrated into a formulation that accounts for both the rigid and flexible behavior of aircraft. A time domain representation of the unsteady aerodynamics is achieved by rational function approximation of the fully unsteady aerodynamics obtained via the doublet lattice method. The case study consists of the so called Body Freedom Flutter research model in two different configurations with one or two tanks partially filled with liquid with a mass comprising 25% of the aircraft structure. The results show that linear sloshing dynamics are able to change the stability margin of the aircraft in addition to having non-negligible effects on rigid body dynamics. On the other hand, vertical sloshing acts as a nonlinear damper and eventually provides limit cycle oscillations after flutter onset.