Maciej Balajewicz


Computational models of high-dimensional systems arise in a rich variety of engineering and scientific contexts. For example, computational fluid dynamics (CFD) and finite element (FE) analysis have become indispensable tools for many engineering applications. Unfortunately, high-fidelity simulations are often computationally prohibitive for parametric and time-critical applications such as design, design optimization, control, and virtual testing. Moreover, even when adequate computational resources are available, simulations often provide too little understanding of the solutions they produce. There are significant scientific and engineering benefits in developing and studying low-dimensional representations of high-dimensional systems that retain physical fidelity while substantially reducing the size and cost of the computational model.

My research is focused on developing theoretical and computational tools for low-dimensional and low-rank models of multi-scale and multi-physics problems.

Data-driven turbulence modeling

  • Balajewicz, M., Tezaur, I., and Dowell, E. “Minimal subspace rotation on the Stiefel manifold for stabilization and enhancement of projection-based reduced order models for the compressible Navier-Stokes equations,” Journal of Computational Physics. Submitted April 2015. (PDF) (arXiv)
  • Balajewicz, M., Dowell, E., and Noack, B. “Low-dimensional modeling of high-Reynolds-number shear flows incorporating constraints from the Navier-Stokes equation,” Journal of Fluid Mechanics. 729:285-308. (2013). (html) (PDF) (arXiv)


Figure 1. Sketch of turbulence by Leonardo da Vinci ca. 1510

It has been said that turbulence is one of the most challenging unsolved problems in Newtonian physics. The aggregate cost to our society resulting from our incomplete understanding of turbulence is significant. Consider, for example, the environmental and economic costs associated with sub-optimal performance of virtually every fluid-thermal system such as the internal combustion engine. It is not possible to directly and routinely simulate turbulent flow because its multi-scale characteristics require prohibitively high computational resources.

My research is an effort to analyze turbulent dynamics in terms of low-dimensional and low-rank direct solutions of the Navier-Stokes equations. To this end, I am the principal author of the successful proposal for a three-year \$479,493 National Science Foundation (NSF) research grant awarded starting August 2014 entitled "Reduced Order Models of the Navier-Stokes Equations of Fluid Flows".

This project focuses on generalizing the proper orthogonal decomposition (POD) to respect the nonlinear energy transfer from productive to dissipative scales, thus creating much more accurate and compact computational models. Conceptually, one might consider this approach as enrichment with additional energy-absorbing dynamic equations, as opposed to empirical eddy-viscosity closure terms. This approach can be expressed mathematically as a constrained optimization problem on the Stiefel manifold.

Preliminary work in this area has been published in the Journal of Fluid Mechanics and also in the journal Nonlinear Dynamics. Future work will involve developing analytical and computational tools for low-dimensional modeling of uncertainty propagation in high-Reynolds-number flows. This focus is motivated by the need to improve fundamental understanding of nonlinear energy transfer between resolved deterministic dynamical components and truncated dynamical components that are characterized by broad energy spectra.

These timely and relevant studies, if successful, will provide a first principles methodology for computing turbulent flows. Firstly, the proposed low-dimensional modeling methodology will have direct and significant impact on numerous practical applications such as active control, bifurcation analysis and optimization that rely heavily on accurate and efficient realizations of turbulent flows. Secondly, the research will shed light on important scientific questions concerning the fundamental nature of turbulence. This will be possible because the low-dimensional models are free of empirical closure models and therefore, constitute direct numerical simulations of turbulence. The low-dimensional models can therefore be used to study the dynamics of the Navier-Stokes equations which are the fundamental equations governing all turbulence. And indeed the methods can ultimately be applied to many other physical systems which have a wide range of spatial and temporal scales in a random like or chaotic response.

Figure 2. Model reduction of high-Reynolds-number flows. (a) Spectral direct numerical simulation (DNS) of two-dimensional lid-driven cavity at Reynolds number at Re = 1 x 105; (b) Temporal evolution of the first modal coefficient and (c) its power spectral density. DNS simulation with 5.2 x 105 DOF (gray lines), standard POD-based model with 10 DOF (dotted line), proposed new model with 10 DOF (red lines)

Machine learning of non-smooth dynamical systems

  • Balajewicz, M., Amsallem, D., and Farhat, C. “Projection-based model reduction for contact problems,” International Journal for Numerical Methods in Engineering. Accepted for publication. (html) (PDF) (arXiv)
  • Balajewicz, M., and Farhat, C. “Reduction of nonlinear embedded boundary models for problems with evolving interfaces,” Journal of Computational Physics. 274:489-504 (2014). (html) (PDF) (arXiv)


Figure 3. Underbody blast

Problems with evolving domains and discontinuities, such as those arising in multi-phase flows, fluid-structure interactions, and contact have inspired the development of a variety of computational methods including embedded or immersed boundary methods, extended finite element methods, and particle methods. Developing model reduction methods for computational frameworks based on such methods seems however to be challenging. Indeed, most popular compression algorithms such as proper orthogonal decomposition (POD) are only applicable to smooth dynamical systems.

This research effort adopts a Netflix-style "big data" approach for tackling reduction of complex non-smooth dynamical systems via clustering, feature recognition, and sparse dictionary learning. In work published in the Journal of Computational Physics, weighted low-rank factorization (WLRA) is utilized to optimally and efficiently model the evolution of sharp discontinuities in highly-nonlinear multi-phase and fluid-structure interaction problems in the embedded boundary CFD framework. In work currently under review for the journal Computer Methods in Applied Mechanics and Engineering, non-negative matrix factorization (NNMF) is adopted to efficiently reduce Karush-Kuhn-Tucker (KKT) multipliers in large-scale constrained optimization problems in contact mechanics. Moreover, KKT-based error indicators are introduced to enable construction of parametrically robust low-dimensional models. In all instances, computational speed-ups of O(103) are reported.

Future work will focus on modeling multi-physics problems with topological changes such as cracking and failure. Projection-based models are global by construction and therefore lack the localization necessary to model inherently local phenomena. For these problems I envision a hybrid approach where local interactions are modeled using compressed sensing and evaluated on a surrogate reduced computational mesh.

Figure 4. Machine learning of non-smooth dynamical systems. (a) Nonlinear finite element (FE) simulation of generic hull underbody blast using 1.4 x 106 DOF. (b) Reduction of model underbody blast problem. High-fidelity FE simulation (black line); standard ROM using singular value decomposition (SVD) (dashed blue lines), proposed ROM using non-negative matrix factorization (NNMF) (dashed red lines). (c) Reduction of model underbody blast problem - self contact between two parallel plates. High-fidelity FE simulation (solid lines), ROM using NNMF (dashed lines). Top floor plate (red lines) bottom floor plate (blue lines)

Fluid-structure interactions using embedded boundary framework

  • Balajewicz, M., and Farhat, C. “Reduction of nonlinear embedded boundary models for problems with evolving interfaces,” Journal of Computational Physics. 274:489-504 (2014). (html) (PDF) (arXiv)


Figure 5. EBM framework

Two common computational frameworks for the solution of fluid-structure interaction problems are the arbitrary Lagrangian Eulerian (ALE) and embedded boundary method (EBM) frameworks. In an ALE framework, the computations are performed on an interface-fitted mesh which is deformed using a mesh motion algorithm to follow the evolution of the interface. In an EBM framework, the interface is embedded in a background mesh - also called an embedding mesh - and allowed to intersect it as it evolves. The EBM framework alleviates many computational challenges, including those associated with meshing complex geometries and is often indispensable for problems where the interface undergoes large-amplitude motions and deformations, and/or topological changes.

Under active development is a fully-implicit Cartesian ghost-fluid EBM solver for the unsteady compressible Navier-Stokes equations. Current capabilities include viscous and inviscid fluid-structure interaction of multiple elastic bodies of arbitrary geometric complexity.

Future work will involve implementation of shared memory multiprocessing capabilities using OpenMP, and feature-driven adaptive mesh refinement.

The main goals of this work are (1) to lower the barrier to entry to numerical research in fluid dynamics, (2) provide short, readable, and easily-modifiable CFD code to enable fast testing of new ideas, (3) serve as a learning platform for undergraduate and graduate students in scientific computing and machine learning.

Figure 6. Fluid-structure interactions using the ghost-fluid embedded boundary method (EBM) framework. (a) Inviscid flow around heaving NACA0012 airfoil at transonic (M = 0.88) flow conditions; Mach number contours. (b) Limit cycle oscillations (LCO) of elastically supported cylindrical body in viscous flow at Reynolds number 500; pressure fluctuation contours. (c) Flapping two-dimensional ellipsoidal body in viscous flow at Reynolds number 1 x 103; vorticity contours

Sparse nonlinear identification of aeroelastic systems

  • Balajewicz, M., and Dowell, E. “Reduced order modeling of flutter and limit-cycle oscillations using the sparse Volterra series,” Journal of Aircraft, 49(6), pp. 1803-1812. (2012). (html) (PDF)
  • Balajewicz, M., Nitzsche, F., and Feszty, D. “Application of multi-input Volterra theory to nonlinear, multi-degree-of-freedom aerodynamic systems,” AIAA Journal, 48(1), pp. 56-62. (2010). (html) (PDF)

The Volterra series has been successfully utilized for the purpose of linearized aeroelastic analysis for over two decades. It has been less successful, however, when applied to nonlinear aeroelastic phenomena such as limit cycle oscillations (LCO) caused by separated flow. The Volterra series - in the aeroelastic context - is a polynomial generalization of Theodorsen's function: the true nonlinear aerodynamic response is approximated progressively using a summation of quadratic, cubic, and higher-order convolution integrals. Unfortunately, since Volterra kernels are multi-dimensional, identification of higher-order kernels using computational and/or experimental training simulations is prohibitive. Thus, even though it is well known that many aerodynamic nonlinearities are poorly approximated by quadratic Volterra models, cubic and higher-order approximations can not be utilized because their identification costs are too high.

In this project we explore novel black-box system identification strategies for highly-nonlinear, multi-input multi-output (MIMO) aeroelastic systems. The key new insight is that for a large class of aeroelastic systems, higher-order Volterra kernels are sparse. Consequently, this sparsity can be exploited to allow for the identification of higher-order approximations at no additional cost. In work published in the AIAA Journal and also in the Journal of Aircraft, 5th order (quintic) Volterra models are demonstrated to successfully predict LCO of the transonic NACA0012 Benchmark Model.

In future work, compressive system identification (CSI) will be adopted to more rigorously identify and optimize sparse realizations of aeroelastic Volterra kernels. Methods developed in this project will be used to model the highly-nonlinear aeroelastic response of flapping-wing micro air vehicles.

Figure 7. Nonlinear system identification of fluid-structure interactions. (a) Inviscid simulation of NACA0012 benchmark model at transonic flow conditions (M=0.8), (b) Limit cycle oscillation (LCO) amplitudes of benchmark model as predicted by high-fidelity CFD (gray lines) and 5th order sparse Volterra series model (dashed lines)