2.1. A note on data and domains

The interdisciplinary nature of PEML can lead to confusion regarding the terms defining the data and domain for the model. To enhance clarity for readers with diverse backgrounds, we provide clarifications on the nomenclature in this paper. Firstly, data refers to all measured or known values that are used in the overall architecture/methodology, not exclusively that which is used as inputs to the ML model, or as target observations. This may include, but is not limited to, measured data, system parameters, and scaling information. In the context of the example above, the data encompasses measured values of the state, $ {\mathbf{z}}^{\ast } $ , and force, $ {f}^{\ast } $ , along with system parameters $ m,c,k,{k}_3 $ (where the asterisk denotes observations of a value). Importantly, the use of the term observation data is akin to the classic ML definition of training data, which is the scope of data used in traditional learning paradigms which minimize the discrepancy between the model output and some observed target values. This change is employed here as the training stages in many instances of PEML demonstrate the learner’s ability to make predictions beyond the scope of these observations of target values.

This goal of extending the scope of prediction also prompts a clarification of the term domain. The domain here is similar to the definition of the domain of a function, representing the set of values passed as input to the model—in this case, the set of time values, $ t $ . The domain where measured values of the model output are available, is termed the observation domain, $ {\Omega}_o $ . The overall domain in which the model is trained, and predictions can be made, is the collocation domain, $ {\Omega}_c $ . For example, if one provides measurements of the state for the first third of the signal in Figure 2b but proposes the model to learn (and therefore predict) over the full signal range, the observation domain would be $ {\Omega}_o\in 0\le t\le 40s $ and the collocation domain $ {\Omega}_c\in 0\le t\le 120s $ . It is crucial to note that these domains are not restricted to a range (scope), and the discrete nature of the observation domain influences the motivation for interpolation schemes. For example, in sparse data recovery schemes, the observation domain can be defined in a discrete manner. Figure 3 provides a visualization of commonly used domain types for a selection of schemes which can employ PEML methods.

Figure 3. Visualization of domain definitions for schemes and motivations that can employ PEML. The blue areas represent the continuous collocation domain, and the red dots represent the coverage and sparsity of the discrete observation domain. The dashed and solid lines represent the scope of the collocation and observation domains, respectively.

5.1. ML-enhanced Bayesian filtering

As previously stated, classical Bayesian filtering requires the model form to be known a priori, implying that the resulting accuracy will depend on how exhaustive this model is. To overcome the inaccuracies that result from model-system discrepancy, ML can be infused with BF techniques to improve inference potential. To this end, Tatsis et al. (Reference Tatsis, Agathos, Chatzi and Dertimanis2022) propose to fuse BF with a Covariance Matrix Adaptation scheme, to extract the unknown position and location of flaws in the inverse problem setting of crack identification, while simultaneously achieving virtual sensing. The latter is the outcome of a hierarchical BF approach powered by reduced order modeling.

Using a different approach, Revach et al. (Reference Revach, Shlezinger, Ni, Escoriza, Van Sloun and Eldar2022) employ a neural network within a Kalman filter scheme to discover the full form of partially known and observed dynamics. By exploiting the nonlinear estimation capabilities of the NN, they managed to overcome the challenges of model constraint, that are common in filtering methods (Aucejo et al., Reference Aucejo, De Smet and Deü2019). Using a similar approach, but with a different motivation, Angeli et al. (Reference Angeli, Desmet and Naets2021) combined Kalman filtering with a deep-learning architecture to perform model-order reduction, by learning the mapping from the full-system coordinates, to a minimal coordinate latent space.

5.2. Physics-guided neural networks

In physics-guided machine learning (PGML), deep learning techniques are employed to capture the discrepancy between an explicitly defined model based on prior knowledge and the true system from which data is attained. The goal is to fine-tune the overall model’s parameters (i.e., the prior and ML model) in a way that the physical prior knowledge steers the training process toward the desired direction. By doing so, the model can be guided to learn latent quantities that align with the known physical principles of the system. This ensures that the resulting model is not only accurate in its predictions but also possesses physically interpretable latent representations. At this stage, we would like to remind the reader of the definitions of physics model strictness and physics constraint given in Section 1. PGML approaches employ a relatively high reliance on the physics-based model form, in that the assumed physics is imposed in a strict form. However, in order to allow for simulation of model discrepancy, the level of physics constraint is relatively low, which means that the learner is not forced to narrowly follow this assumed prescribed form. The relaxation of constraints differs between PGML and physics-informed (PIML) approaches; PGML techniques often reduce constraints through the use of bias or residual modeling, whereas PIML schemes employ physics into the loss function as a target solution (i.e., the solutions are weakly imposed).

One of the key advantages of physics-guided machine learning is its ability to incorporate domain knowledge into the learning process. This is particularly beneficial in scenarios where data may be limited, noisy, or expensive to obtain. In comparison to PEML techniques that lay lower on our prescribed spectrum, PGML methods can offer an increased level of interpretability. By steering the model’s learning process with physics-based insights, it can more effectively generalize to unseen data and maintain a coherent understanding of the underlying physical mechanisms.

The training process in physics-guided models involves two key components: (1) Incorporating Prior Knowledge: Prior knowledge on the physics of the system is integrated into the network architecture, or as part of the model; (2) Capturing Discrepancy: Deep learning models excel in learning from data, even when this contradicts prior knowledge. As a result, during training, such a model will gradually adapt and learn to account for discrepancies between the prior knowledge and the true dynamics of the system. This adaptability allows the model to converge toward a more accurate representation of the underlying physics.

Conceptually similar to estimating a residual modeler, Liu et al. (Reference Liu, Lai and Chatzi2021, Reference Liu, Lai, Bacsa and Chatzi2022) proposed a probabilistic physics-guided framework termed a Physics-guided Deep Markov Model (PgDMM) for inferring the characteristics and latent structure of nonlinear dynamical systems from measurement data. It addresses the shortcoming of black-box deep generative models (such as the DMM) in terms of lacking physical interpretation and failing to recover a structured representation of the learned latent space. To overcome this, the framework combines physics-based models of the partially known physics with a DMM scheme, resulting in a hybrid modeling approach. The proposed framework leverages the expressive power of DL while imposing physics-driven restrictions on the latent space, through structured transition and emission functions, to enhance performance and generalize predictive capabilities. The authors demonstrate the benefits of this fusion through improved performance on simulation examples and experimental case studies of nonlinear systems.

Both residual modeling and PgNNs, in general, share a common objective of easing the training objective of neural networks. Residual modeling achieves this by easing the learning process through a general approximation, while PgNNs incorporate physical knowledge to provide reliable predictions even in data-limited situations. However, regarding specific implementations, there are differences between these two methods. Residual modeling operates without the need for specific domain knowledge, focusing instead on the residuals in a more general-purpose application. This makes it broadly applicable across various standard machine learning tasks, especially in the domain of computer vision. PgNNs, on the other hand, explicitly incorporate physical laws into the model, adding interpretability related to these physical models. PgNNs can also work in a general setting with any prior models to fit the residuals, but the key idea is to use physical prior model to obtain a physically interpretable latent representation from neural networks. This makes PgNNs particularly suited for problems where adherence to physical laws is paramount. Additionally, while residual networks rely heavily on data for learning, PgNNs can leverage physical models to make predictions even with limited data, demonstrating their utility in data-constrained environments.

To demonstrate how a physical prior model can impact the training of DL models, we apply the PgDMM to our working example, using the same data generation settings as explained in Section 4. In this case, instead of approaching the system with no prior information, we introduce a physical prior model into the DMM to guide the training process. This physical prior model is a linear model that excludes the cubic term in equation (1), which replicates a knowledge gap in the form of additional system complexities. The results are shown in Figure 6. It can be observed that the predictions for both displacements and velocities align well with the ground truth. The estimation uncertainty is slightly higher for velocities, which is expected since they are unobserved quantities. It is important to note that the system displays significant nonlinearity due to the presence of a cubic term with a large coefficient, causing the linear approximation to deviate noticeably from the true system dynamics. However, the learning-based model within the framework still captures this discrepancy and reconstructs the underlying dynamics through the guided training process.

Figure 6. Predictions versus exact solutions of displacement (top) and velocity (bottom) using the PgDMM applied to the working example. Displacement is assume to be the only measurement. The gray dash-dot line is the physical prior model and the blue bounding boxes represent the estimated $ 2\sigma $ range.

A similar physics-guided RNN was proposed by Yu et al. (Reference Yu, Yao and Liu2020), which consists of two parts: physics-based layers and data-driven layers, where physics-based layers encode the underlying physics into the network and the residual block computes a residual value which reflects the consistency of the prediction results with the known physics and needs to be optimized toward zero.

Instead of modeling the residual of the prior model, a physics-guided Neural Network (PgNN), proposed by Karpatne et al. (Reference Karpatne, Watkins, Read and Kumar2017), not only ingests the output of a physics-based model in the neural network framework, but also uses a novel physics-based learning objective to ensure the learning of physically consistent predictions, as based on domain knowledge. Similarly, the authors proposed a Physics-guided Recurrent Neural Network scheme (PGRNN, Jia et al., Reference Jia, Willard, Karpatne, Read, Zwart, Steinbach and Kumar2019) that contains two parallel recurrent structures—a standard RNN flow and an energy flow to be able to capture the variation of energy balance over time. While the standard RNN flow models the temporal dependencies that better fit observed data, the energy flow aims to regularize the temporal progression of the model in a physically consistent fashion. Furthermore, in another PgNN proposed by Robinson et al. (Reference Robinson, Pawar, Rasheed and San2022), the information from the known part of the system is injected into an intermediate layer of the neural network.

The physics-guided deep neural network (PGDNN), proposed by Huang et al. (Reference Huang, Yin and Liu2022), uses a cross-physics-data domain loss function to fuse features extracted from both the physical domain and data domain, which evaluates the discrepancy between the output of a FE model and the measured signals from the real structure. With the physical guidance of the FE model, the learned PGDNN model can be well generalized to identify test data of unknown damages. The authors also use the same idea in bridge damage identification under moving vehicle loads (Yin et al., Reference Yin, Huang and Liu2023). Similarly, Zhang and Sun (Reference Zhang and Sun2021) presented usage of the FE model as an implicit representation of scientific knowledge underlying the monitored structure and incorporates the output of FE model updating into the NN model setup and learning.

In Chen and Liu (Reference Chen and Liu2021), the physics knowledge is incorporated into the neural network by means of imposing appropriate constraints on weights, biases or both. Muralidhar et al. (Reference Muralidhar, Bu, Cao, He, Ramakrishnan, Tafti and Karpatne2020) used physics-based prior model, physical intermediate variables, and physics-guided loss functions to learn physically interpretable quantities such as pressure field and velocity field.

6.1. Dictionary methods

One of the biggest challenges faced in the practical application of structural mechanics in engineering, is the presence of irregular, unknown, or ill-defined nonlinearities in the system. Another challenge may occur from variation in the parameters which govern the prescribed model of the system, which could be from environmental changes, or from consequential changes such as damage. This motivates less reliance on the physics-based model form, to allow for freedom in physical-digital system discrepancy whilst satisfying known physics. One approach to reducing the reliance on the prescribed physics model form is by having the learner estimate the definition of the model, which may be the sole, or additional, objective of the learner. Dictionary methods are well-positioned as a solution for less strict physics embedding, where the model is determined from a set of possible model solutions, allowing freedom in a semi-discrete manner.

The problem of estimating the existence, type, or strength of the model-governing physics is described as that of equation discovery. This inverse problem can often be very computationally expensive, due to the large number of forward model calculations required to evaluate the current estimate of the parameters (Frangos et al., Reference Frangos, Marzouk, Willcox and van Bloemen Waanders2010). When determining the presence of governing equation components, often the identification is drawn from a family of estimated equations. For example, the matching-pursuit algorithm selects the most sparse representation of a signal from a dictionary of physics-based functions (Vincent and Bengio, Reference Vincent and Bengio2002).

For dictionary-based approaches, the idea is to determine an estimate of the model output as some combination of bases or “atoms.” Often, these bases are formed as candidate functions of the input data, and are compiled into a dictionary-matrix $ \Theta \left(\mathbf{x}\right) $ . Then, linear algebra is used to represent the target signal from a sparse representation of this dictionary and a coefficient matrix $ \Xi $ .

(8) $$ \mathbf{y}=\Theta \left(\mathbf{x}\right)\Xi . $$

And so, the aim of the learner is to determine a suitably sparse solution of $ \Xi $ , via some objective function. Typically, sparsity-promoting optimisation methods are used, such as LASSO (Tibshirani, Reference Tibshirani1996). The coefficient matrix can be constrained to be binary values, thus operating as a simple mask of the candidate function, or can be allowed to contain continuous values, and thus can simultaneously determine estimated system parameters which are used in the candidate functions. For the case of dynamic systems, a specific algorithm was developed by Brunton et al. (Reference Brunton, Proctor and Kutz2016a) for sparse discovery of nonlinear systems. The team also showed how this could be used in control by including the control parameters in the dictionary definition (Brunton et al., Reference Brunton, Proctor and Kutz2016b). Kaiser et al. (Reference Kaiser, Kutz and Brunton2018) extended this nonlinear dictionary learning approach to improve control of a dynamic system where data is lacking. To do so, they extended the method to include the effects of actuation for better forward prediction.

In an example of dictionary-based learning, Flaschel et al. (Reference Flaschel, Kumar and De Lorenzis2021) and Thakolkaran et al. (Reference Thakolkaran, Joshi, Zheng, Flaschel, De Lorenzis and Kumar2022) showed a method for unsupervised learning of the constitutive laws governing an isotropic or anisotropic plate. The approach is not only unsupervised, needing only displacement and force data, it is directly inferrable in a physical manner. The authoring team then extended this work further to include a Bayesian estimation (Joshi et al., Reference Joshi, Thakolkaran, Zheng, Escande, Flaschel, De Lorenzis and Kumar2022), allowing for quantified uncertainty in the model of the constitutive laws.

Another practical example: Ren et al. (Reference Ren, Han, Yu, Skjetne, Leira, Sævik and Zhu2023) used nonlinear dynamic identification to successfully predict the forward behavior of a 6DOF ship model, including coupling effects between the rigid body and water. In this work, the dictionary method is combined with a numerical method to predict the state of the system in a short time window ahead. This facet is often found when applying DM-based approaches to practical examples, as the method is intrinsically a model discovery approach, and so a solution step is required if model output prediction is required.

Data-driven approaches to equation discovery will often require an assumption of the physics models being exhaustive of the “true” solution; that is all the parameters being estimated will fully define the model, or the solution will lie within the proposed family of equations. One of the challenges faced with deterministic methods, such as LASSO (Tibshirani, Reference Tibshirani1996), is their sensitivity to hyperparameters (Brunton et al., Reference Brunton, Proctor and Kutz2016a); a potential manifestation of which is the estimation of a combination of two similar models, which is a less accurate estimate of the “true” solution than each of these models individually. Bayesian approaches can help to overcome this issue, by instead providing a stochastic estimate of the model, and enforcing sparsity (Park and Casella, Reference Park and Casella2008).

Fuentes et al. (Reference Fuentes, Nayek, Gardner, Dervilis, Rogers, Worden and Cross2021) show a Bayesian approach for nonlinear dynamic system identification which simultaneously selects the model, and estimates the parameters of the model. Similarly, Nayek et al. (Reference Nayek, Fuentes, Worden and Cross2021) identify types and strengths of nonlinearities by utilizing spike-and-slab priors in the identification scheme. As these priors are analytically intractable, this allowed them to be used along with a MCMC sampling procedure to generate posterior distributions over the parameters. Abdessalem et al. (Reference Abdessalem, Dervilis, Wagg and Worden2018) showed a method for approximate Bayesian computation of model selection and parameter estimation of dynamic structures, for cases where the likelihood is either intractable or cannot be approached in closed form.

6.2. Physics-informed neural networks

Raissi et al. (Reference Raissi, Perdikaris and Karniadakis2019) showed that by exploiting automatic differentiation that is common in the practical implementation of neural networks, one can embed physics that is known in the form of ordinary/partial differential equations (ODEs/PDEs). Given a system of known ODEs/PDEs which define the physics, where the sum of these ODEs/PDEs should equal zero, an objective function is formed which can be estimated using automatic differentiation over the network. If one was to apply a PINN to estimate the state of the example in Figure 2, over the collocation domain $ {\Omega}_c $ , using equation (2), the physics-informed loss function becomes,

(9) $$ {L}_p\left(\mathbf{t};\Theta \right)={\left\langle {\partial}_t{\mathcal{N}}_{\mathbf{z}}-\mathbf{A}{\mathcal{N}}_{\mathbf{z}}-{\mathbf{A}}_n{\mathcal{N}}_{u^3}-\mathbf{B}f\right\rangle}_{\Omega_c}, $$

where $ {\partial}_t{\mathcal{N}}_{\mathbf{z}} $ is the estimated first-order derivative of the state using automatic differentiation, and the system parameters are $ \theta =\left\{m,c,k,{k}_3\right\} $ . Boundary conditions are embedded into PINNs in one of two ways; the first is to embed them in a “soft” manner, where an additional loss term is included based on a defined boundary condition $ \xi \left(\mathbf{y}\right) $ (Sun et al., Reference Sun, Gao, Pan and Wang2020). Given the boundary domain $ \mathrm{\partial \Omega}\in {\Omega}_c $ ,

(10) $$ {L}_{bc}={\left\langle {\xi}^{\ast }-\xi \left({\mathcal{N}}_{\mathbf{y}}\right)\right\rangle}_{\mathrm{\partial \Omega }}. $$

The second case, so-called “hard” boundary conditions, involves directly masking the outputs of the network with the known boundary conditions,

(11) $$ {\left.{\mathcal{N}}_{\mathbf{y}}={\xi}^{\ast}\right|}_{\mathrm{\partial \Omega }}, $$

where often gradated masks are used to avoid asymptotic gradients in the optimization (Sun et al., Reference Sun, Gao, Pan and Wang2020). In the case of estimating the state over a specified time window, the boundary condition becomes the initial condition (i.e., the state at $ t=0 $ ).

(12) $$ \xi =\mathbf{z}(0)={\left\{u(0),\dot{u}(0)\right\}}^T,\hskip2em {L}_{bc}={\left\langle \xi -{\mathcal{N}}_{\mathbf{z}}\right\rangle}_{\Omega \in t=0}. $$

The total objective function of the PINN is then formed as the weighted sum of the observation, physics, and boundary condition losses,

(13) $$ L={\lambda}_o{L}_o+{\lambda}_p{L}_p+{\lambda}_{bc}{L}_{bc}. $$

The flexibility of the PINN can then be highlighted by considering how the belief of the architecture is changed when selecting the corresponding objective weighting parameters $ {\lambda}_k $ . By considering the weighting parameters as 1 or 0, Table 1 shows how the selection of objective terms to include changes the application type of the PINN. The column titled $ {\Omega}_c $ indicates the domain used for the PDE objective term, where the significant difference is when the observation domain is used, and therefore the PINN becomes akin to a system identification problem solution.

Table 1. Summary of PINN application types, and the physics-enhanced machine learning genre/category each would be grouped into

Figure 7 shows the framework of a generic PINN, highlighting where the framework, specifically the loss function formulation, can be broken down into the data-driven and physics-embedding components. It is also possible to include the system parameters $ \theta $ as unknowns which are determined as part of the optimisation process, that is, $ \Theta =\left\{\mathbf{W},\mathbf{B},\theta \right\} $ . Doing so adds the capability as an system identification tool, as well as equation solution discovery. In system identification, either equation discovery can be performed, where the aim is to determine the definition of $ \mathrm{\mathcal{F}} $ , or parameter estimation, where the aim is to determine the values of $ \theta $ . The framework shown is for the training stage of the model, where the network parameters $ \Theta $ , and optionally the physical parameters $ \theta $ , are updated using an optimization algorithm such as LBFGS (Liu and Nocedal, Reference Liu and Nocedal1989) or Adam (Kingma and Ba, Reference Kingma and Ba2014). Further prior physical knowledge can be added to PINNs by embedding more objective functions, such as initial conditions or continuity conditions, which can be formed directly from the output values of the network, or from the derivatives.

Figure 7. Framework of a general PINN, highlighting where the data-driven and physics-knowledge are embedded within the process.

To illustrate the different approaches to implementing a PINN, three paradigms of the method are discussed here, each to aiming to tackle a different objective.

1. 1. System Identification; the aim is to determine the physics description of the system, either by equation discovery or by parameter estimation of the system parameters $ \theta $ in equation (9) via estimating the system state.

2. 2. Enhanced Learning; the aim is to enhance the learner to either better estimate the model given a more sparse set of data, or to improve learning efficiency.

3. 3. Forward Modeling; where the PINN acts to generate a “simulated” model of the system, given the system equations and parameters, and the domain of interest.

These points will form the remaining parts of this subsection, where surveys will be shown and an example of each approach applied to the working example of the Duffing oscillator. The PINN is applied here as an “instance-modeler”; where the estimated model is only applicable in the case of the prescribed initial conditions and forcing signal. Therefore, new estimations made with this model are only applicable within the training domain. The generalisability and extrapolability presented with the PINN are in terms of extending beyond the domain of observations.

Physics-informed neural networks are gaining increased attention thanks to some of the advantages stated above, however, there exist a number of drawbacks which should be noted. Firstly, in many formulations of PINNs, they suffer from a lack of generalisability—which is a common motivation for PEML schemes—as they are restricted to the domain on which they were trained (Haghighat et al., Reference Haghighat, Raissi, Moure, Gomez and Juanes2021b). This domain may be extended beyond that of observations, however, computationally intensive training must still be performed before prediction. Furthermore, there is a lack of intuition or knowledge on the optimization task of PINNs; often, weighting of the losses is done via trial and improvement (Wang et al., Reference Wang, Teng and Perdikaris2021), and an under-defined, or ill-posed, physics prescription may result in many local minima (Nandi et al., Reference Nandi, Hennigh, Nabian, Liu, Woo, Jordan, Shahnam, Syamlal, Guenther and VanEssendelft2021). A related challenge is in the computational effort of training, which can be a considerable task. This challenge has garnered fair criticism against the use of PINNs as opposed to other numerical-solving schemes (Grossmann et al., Reference Grossmann, Komorowska, Latz and Schönlieb2023).

6.2.1. PINNs for system identification

As discussed above, PINNs can be used as an approach to system identification, by applying a “soft” condition on the governing physical laws. The use of a variety of soft conditions, in the form of the loss functions, allows for discrepancies between the model and data, making them useful for noisy observations. A good demonstration of PINNs for simultaneous state-system estimation is shown by Yuan et al. (Reference Yuan, Zargar, Chen and Wang2020) and Moradi et al. (Reference Moradi, Duran, Eftekhar Azam and Mofid2023), where in both examples, they model the displacement of a vibrating beam, with an accurate estimation of the governing equation parameters also. Zhang et al. (Reference Zhang, Yin and Karniadakis2020) use a PINN to solve an identification problem for nonhomogeneous materials, using elasticity imaging. An elegant approach here was done so by utilizing a multi-network architecture to include nonhomogeneous parameter fields, removing the potential issue caused by exploding dimensionality when including the spatially dependent material parameters. Sun et al. (Reference Sun, Liu, Wang and Sun2023) recently showed an approach to discovering the parameters of complex nonlinear dynamic systems using a Physics-informed Spline Learning approach. This approach employs the same exploitation as PINNs, however, the splines are used to allow for differentiation from a more sparse set of observed data, by interpolating the underlying dynamics.

Instead of determining the governing parameters of the PDEs, Wu and Xiu (Reference Wu and Xiu2020) proposed a method to estimate the forward-in-time solution to a model. By recovering the evolution operator in the modal space, this reduces the problem from an infinite-dimension to a finite-dimension space. As opposed to value-based parameter identification, Ritto and Rochinha (Reference Ritto and Rochinha2021) used measurement from a vibrating bar to update a ML-based classifier which directly infers the damage state of the structure. As the governing equations are often in the form of PDEs, an estimation of the derivatives is required in the learning process of many equation discovery approaches. However, Goyal and Benner (Reference Goyal and Benner2022) proposed a numerical integration framework with dictionary learning, along with a “Runge–Kutta inspired” numerical procedure, overcoming the issues presented with derivative approximation from corrupted or sparse data.

For the working example, the PINN was applied to the Duffing oscillator, shown in Figure 2, and an equivalent linear system (i.e., $ {k}_3=0 $ ). Here, 256 data points were sub-sampled using a Sobol sampler and passed to the learner, and the physics-loss domain is set to the same as the observation domain (i.e., $ {\Omega}_p={\Omega}_o $ ), positing the framework akin to a system identification scheme. The observation loss (equation (7)) penalizes the output enough to result in an accurate estimation of the displacement, and satisfying the physics loss will drive the estimation of the physical parameters.

The results of the state estimation are shown in Figure 8, and the results of the parameter estimation are shown in Table 2. The state estimation results show accurate modeling, and the estimated values for the physical parameters also show a good level of accuracy. It is important to note here that accurate estimation of the state is not the primary objective, given that the domain of interest is well-covered by the observation data. A notable result is the increased accuracy when modeling only the linear system. This result can be explained by well-known machine learning intuition that with an increased dimensionality of estimation-space, an increased number of information is also required. Therefore, with the same level of information provided, the accuracy of the results for estimating more physical parameters will likely be more of a challenge for the learner.

Figure 8. Predicted versus exact solution of simultaneous system-state estimation approach to solving the working example for the nonlinear case (top) and linear case ( $ {k}_3=0 $ ) (bottom).

Table 2. Results of system estimation for the SDOF oscillator for both the nonlinear and linear case

Note. For all parameters, the estimated value and the percentage error are shown.

6.2.2. PINNs for enhanced learning

Another utility of PINNs is domain enhancement potential, either by improving domain density from sparse observations, or by extending the domain beyond that of the observation domain. Practically for spatio-temporal models, increasing density is often only motivated in the spatial domain, as a result of the ease of improving time-domain sampling. However, domain-extension approaches are found to be practically motivated in both space and time.

Xu et al. (Reference Xu, Han, Cheng, Cheng and Ge2022) used PINNs to accurately model the rigid-body dynamics of an unmanned surface vehicle as it voyages along a river. This is a nice example of a practical implementation of improving state prediction from sparse data, and a good example of how to formulate PINNs for relatively complex descriptions of dynamics. Chen and Liu (Reference Chen and Liu2021) use a PINN for estimating the fatigue $ S-N $ curves, where even on seemingly sufficient data, the uninformed ANN fails to accurately predict. A particular note of this work, was the inclusion of a probabilistic framework, allowing both freedom in the model construction, as a result of stochastic considerations, and a quantified estimate of the uncertainty. By utilizing a finite-element model, in which the parameters were updated using a PINN, Zhang and Sun (Reference Zhang and Sun2021) developed a NN-based method for detecting damage.

By utilizing an energy-based formulation of the loss function, Zhuang et al. (Reference Zhuang, Guo, Alajlan, Zhu and Rabczuk2021) modeled bending, vibration, and buckling of a Kirchoff plate. As well as the informed loss function, the authors used a non-standard activation function to better emulate the underlying operations which govern the physics. Another example of an energy-based loss function approach is shown by Goswami et al. (Reference Goswami, Anitescu, Chakraborty and Rabczuk2020), who combined a PINN with transfer learning to model the phase field of fracture in a material. They showed how a well-trained model could drastically reduce the computational requirements of this problem.

The deformation of elastic plates with PINNs was shown by Li et al. (Reference Li, Bazant and Zhu2021), where they made a comparison between purely data-driven, PDE-based, and energy-based physics informing, finding each to have different advantages. The PDE-based approach was less dependent on sampling size and resolution, whereas the energy-based approach had less hyperparameter, and therefore was more efficient and easier to train. This information may be useful for the future development of PINNs, which suffer from a lack of understanding of the hyperparameter space pathology, to allow for robust optimization strategies, as discussed by Wang et al. (Reference Wang, Teng and Perdikaris2021).

Yucesan and Viana (Reference Yucesan and Viana2020) present a methodology for predicting the damage level in a wind turbine blade, in the form of grease degradation, using PINNs. In their work, the advantage of PINNs is that instead of aiming to directly model the value of grease degradation, they aim to predict the increment of grease degradation based on the current value and a number of other measurable quantities. This would be a difficult task to perform in a black-box manner, as the values of the increment are not the observed values, and so the physics embedding helps to overcome this.

Now, we will show the PINN applied as an enhanced learner to the working example, where the physical parameters are known a priori. The advantages of the PINN in this aspect are more efficient learning, and for improved extrapolation of the data from sparse observations. To demonstrate, only every sixteenth sample is instead fed to the learner, emulating a sampling frequency of 0.5328 Hz. Then, the PINN is applied as a “black-box” architecture (i.e., $ {\lambda}_{ic}={\lambda}_p=0 $ ), and as a physics-informed modeler, where all loss weight parameters are included. The results of the two methods applied are shown in Figure 9, where it is clear to see the strength of embedding physics, demonstrated in the context of sparse data.

Figure 9. Predicted versus exact solution of state estimation approach applied to a subsample of the working example with no physics embedded (top) and physics-informed embedding (bottom).

6.2.3. PINNs for forward modeling

When applying PINNs for forward modeling of systems and structures, no prior observations are given, that is, , however, it is necessary to provide sufficient boundary and initial conditions, and physics constraints which describe a complete model. An advantage of PINNs for forward modeling is their simple implementation; embedding boundary conditions, complex geometries, or new governing equations is relatively straightforward. The training time of PINNs as forward modelers, in relation to traditional finite-element methods, is often greater, leading to statements implying their impracticality (Rezaei et al., Reference Rezaei, Harandi, Moeineddin, Xu and Reese2022). In the context of forward modeling of a known model, and at pre-set collocation points, this statement of impracticality is well-founded. However, PINNs may also provide a convenient solution for problems such as; in-time control/prediction, where a PINN could predict an output rapidly, as the computational effort is done a priori during training, or for determining a better generalization over high-dimensional modeling space, where the computational cost of FEM methods can rapidly increase with dimensionality.

PINNs for forward modeling have gained a lot of traction in micro-scale problems; Haghighat et al. (Reference Haghighat, Raissi, Moure, Gomez and Juanes2021b) and Henkes et al. (Reference Henkes, Wessels and Mahnken2022) use PINNs to model the instance of displacements and stresses in a unit cell. In the latter, they showed the capability of the approach to model nonlinear stresses by including a sharp phase transition within the material. Haghighat and the team also showed how the method could accurately model structural vibration (Haghighat et al., Reference Haghighat, Bekar, Madenci and Juanes2021a), giving only initial and boundary conditions. Abueidda et al. (Reference Abueidda, Lu and Koric2021) showed PINNs for modeling various solid mechanic effects; elasticity, hyperelasticity, and plasticity, where the method performed well on all types of materials. In the work by Zheng et al. (Reference Zheng, Li, Qi, Gao, Liu and Yuan2022), fracture mechanics were modeled to a decent accuracy with only principle physics.

The above examples all follow a fairly straightforward path to the model formulation, by employing the variables that are inherent to the governing equations as the variables of the ML model. Going beyond this, Huang et al. (Reference Huang, Fuhg, Weißenfels and Wriggers2020) used a proper orthogonal decomposition (POD) neural network to model plasticity in a unit cell. The advantage of the POD approach is to decouple the multi-dimensional stress, allowing the use of individual NNs for each stress variable, reducing computation time, and increasing learning efficiency. Straying away from the common approach of using PDEs to form the physics-based loss, Abueidda et al. (Reference Abueidda, Koric, Al-Rub, Parrott, James and Sobh2022) formulated an energy-based loss term to successfully forward model hyperelasticity and viscoelasticity in a given material.

As a final demonstrator of PINN utility, we apply it as a forward modeler, with no observed data offered to the learner, that is, . However, it is important to note that there is still training data reflecting time domain information. For this implementation, initial conditions, in the form of Dirichlet and Neumann boundary conditions from equation (12), and the forcing signal are provided. In Figure 7, therefore, only the physics embedding portion of the framework is included. The results of this forward modeling approach are shown in Figure 10, where it can be seen the solution matches well with the exact solution.

Figure 10. Exact solution versus PINN-based forward modeling solutions of the SDOF Duffing oscillator example, where no observations of the state are given to the learner.

7.1. Constrained GPs

One example of physics-encoded learners is constrained Gaussian processes (GPs), which, depending on kernel design, can be viewed as embedding the general shape of the function, as GPs are a problem of discovery over the function-space, as opposed to the weight-space (O’Hagan, Reference O’Hagan1978; Williams and Rasmussen, Reference Williams and Rasmussen2006). Conceptually, this process can be thought of as estimating a distribution over all the possible functions that could explain the data, as opposed to one “best fit” model. The aim is to estimate a nonlinear regression model given a set of observed output data $ \mathbf{y} $ , and observed input data $ \mathbf{x} $ ,

(14) $$ \mathbf{y}=f\left(\mathbf{x}\right)+\varepsilon, \hskip2em \varepsilon :\mathcal{N}\left(0,{\sigma}_n^2\unicode{x1D540}\right), $$

where $ \varepsilon $ is a zero-mean Gaussian white-noise process with variance $ {\sigma}_n^2 $ . A GP is fully defined by its mean and covariance functions,

(15) $$ f\left(\mathbf{x}\right):\mathcal{GP}\left(m\left(\mathbf{x}\right),K\left(\mathbf{x},{\mathbf{x}}^{\prime}\right)\right). $$

The mean function $ m\left(\mathbf{x}\right) $ can be any parametric mapping of $ \mathbf{x} $ , and the covariance function expresses the similarity between two input vectors $ \mathbf{x} $ and $ {\mathbf{x}}^{\prime } $ . The primary influence of the user, when implementing a GP, is in the choice of the covariance kernel, which is calculated as any other kernel; linear pair-wise distances between points to form a covariance matrix. There are a number of popular kernels that are used as standard, each of which embed a different belief as to which family of functions the model solution is drawn from. There are two primary methods of embedding physical knowledge in GPs; the first of which is to include an initial estimate of the model into the mean function, and so the problem can then be envisioned as determining the solution of the remaining, unknown physics.

The second approach is to design, or select, the covariance kernels to constrain the shape of the function estimate, and often combinations of kernels can provide varying levels of physical knowledge embedding (Cross and Rogers, Reference Cross and Rogers2021). In this case, it is possible to relax the reliance on the physics model form, by for example simply dictating the design of the kernel for the domain on which the physics is expected to operate. For example, Padonou and Roustant (Reference Padonou and Roustant2016) showed a kernel for predicting on circular domains. More strictly prescribed model forms can be applied by designing kernels to include physical knowledge in the form of partial differential equations, or boundary conditions (Solin and Särkkä, Reference Solin and Särkkä2020).

An example of embedding prior knowledge via the mean function is shown by Zhang et al. (Reference Zhang, Rogers and Cross2021), where they showed an improved GP model for the deflection of the Tamar bridge (Cross et al., Reference Cross, Koo, Brownjohn and Worden2013) under time-varying environmental conditions, by including an expected linear deflection of the cable due to temperature. Data from the Tamar bridge represents a relatively complex problem; a large-scale complicated structure under varying (and non-exhaustive) environmental conditions. However, even by applying the simplest method of embedding physics into a GP, the modeling was improved. Petersen et al. (Reference Petersen, Øiseth and Lourens2022) also applied a novel physics-informed GP method to a bridge problem, but with the aim of estimating wind load from acceleration data. In their work, they developed a novel infusion of GP latent-force model (GP-LFM) with a Kalman filter-based approach. The inclusion of the GP-LFM allowed for characterization of the evolution of the wind-load, and this is enriched with prior physical knowledge in the form of stochastic information on wind-loads taken from wind-tunnel tests. This work provides an excellent demonstration of how physics information can be embedded to allow the transfer of information from scaled structures.

Haywood-Alexander et al. (Reference Haywood-Alexander, Dervilis, Worden, Cross, Mills and Rogers2021) used constrained GPs to model the physical characteristics of guided-waves in a complex material. Guided waves in complex materials are famously difficult to model due to their relatively short wavelength in comparison to the material structure. By designing a variety of kernels, they demonstrated the performance of the GP modeling with varying levels of physics information embedded. Notably, they showed that by considering the space in which the physics operates, one can already improve the modeling capabilities, even before any physics equations are embedded. Continuing on the topic of elastic waves, Jones et al. (Reference Jones, Rogers and Cross2023) applied PI-GPs to localize acoustic source emission in a complex domain. They developed constrained GPs further by embedding boundary conditions and the spatial domain of the problem. This approach has potential for use in modeling on structures with relatively high geometrical complexity, such as those with lots of joints, or with layered materials.

As stated above, the GP is unique in that it operates in the function space, thus, prior knowledge embedded is on that of the shape of the function. This characteristic was utilized by Dardeno et al. (Reference Dardeno, Haywood-Alexander, Mills, Bull, Dervilis and Worden2021), who used weak-form dynamics equations as a mean function within a novel overlapping mixture of Gaussian processes (OMGP) method. By constraining the expected shape of the functions, this allows the learner to separate out unsorted data of dynamic structures from within a population.

To demonstrate the constrained GP, here were apply this to the working example in Figure 2. In this context, the GP estimates a nonlinear operator where the input is time (i.e., $ \mathbf{x}=\mathbf{t} $ ). The first kernel we will select is the scaled squared-exponential kernel, where $ l $ is a lengthscale hyperparameter, and $ \alpha $ is the scaling hyperparameter,

(16) $$ {K}_{SE}\left(\mathbf{t},{\mathbf{t}}^{\ast}\right)=\exp \left(-\frac{1}{2{l}^2}{\left(\mathbf{t}-{\mathbf{t}}^{\ast}\right)}^T\left(\mathbf{t}-{\mathbf{t}}^{\ast}\right)\right), $$

which embeds only a belief that the function is smooth. Then, following the work of Cross and Rogers (Reference Cross and Rogers2021), if we assume a Gaussian white noise force input, an additional physics-derived kernel can be included,

(17) $$ {K}_{\mathrm{SDOF}}\left(\mathbf{t},{\mathbf{t}}^{\ast}\right)=\frac{\sigma_f^2}{4{m}^2{\zeta \omega}_n^3}{e}^{-{\zeta \omega}_n\mid \tau \mid}\left(\cos \left({\omega}_d\tau \right)+\frac{{\zeta \omega}_n}{\omega_d}\sin \Big({\omega}_d|\tau |\Big)\right),\hskip1em \tau =\mathbf{t}-{\mathbf{t}}^{\ast }, $$

where $ {\sigma}_f $ is an additional hyperparameter, representing the square root of the variance of the forcing. To demonstrate the constrained GP approach, the data from Figure 2 is applied in a similar fashion to the informed-learner approach used with the PINN in Figure 9. Here, only every twelfth sample is fed to the learner, and is first applied with the uninformed squared-exponential kernel in equation (16), and then with the constrained SDOF kernel in equation (17). The results of the estimated solutions using each of these kernels are shown in Figure 11, along with the estimated 95% ( $ 2\sigma $ ) confidence range. Note, the value of $ \sigma $ here is taken from the estimated covariance matrix, and is not the hyperparameter $ {\sigma}_f $ . It is clear to see the improved estimation of the displacement of the system, but on top of this, the estimated uncertainty is also reduced. An interesting observations from this is that when using a linear SDOF kernel to model a nonlinear system, this still results in improved modeling from sparse data.

Figure 11. Predicted versus exact solutions of displacement estimation using a GP applied to a subsample of the working example, with (top) no physics embedded and (bottom) constrained GP. The blue bounding boxes represent the estimated $ 2\sigma $ range.

Constrained Gaussian processes are often used to improve the forecasting of temporal data, which falls under the umbrella of domain extension schemes. In the context of the example prescribed here, for accurate forecasting one would need to determine an analytical solution for the covariance of a forced nonlinear system, as with the prescribed kernel, prediction beyond the scope of the data results in simply a free-vibration system.

7.2. Physics-encoded neural networks

So far we have only explored how the automatic differentiation mechanisms can help embed physics into a neural networks. This process can be applied blindly to any type of neural network architecture. A different approach is to modify the architecture of the neural network itself such that its hidden states conform to a domain closer to that of the physical problem of study. By imposing specific geometric constraints, it is possible to bias neural networks toward such domains. Adding additional biases to the model reduces the variance during the training phase, which leads to a faster convergence. However, such biases must be carefully imposed such that they replicate the physics of the system of study, else they add a constant bias error into the model.

The biases that one can introduce in a machine learning model are in most cases symmetries that are imposed on the nonlinear mappings of the model (Bronstein et al., Reference Bronstein, Bruna, Cohen and Veličković2021). The rational for imposing symmetries stems from Noether’s theorem: every continuous symmetry in a physical system corresponds to a conserved quantity. Therefore, symmetries in a model’s architecture should be able to encode the physical properties of the system that are conserved. From henceforth, we will refer to models as “physics-encoded” whenever their architecture is biased to specifically reproduce the symmetries arising from the properties of the system of study. Such symmetries can be tailored depending on the prior knowledge available on the system.

7.2.2. Hamiltonian neural networks

An alternative method for embedding physics into the architecture of the network is to constraint them according the Hamiltonian formalism. The intuition behind this method originates from Noether’s first theorem which states that every differentiable symmetry of the action of a physical system with conservative forces has a corresponding conservation law. In layman terms, this means that the symmetries that are observed within the dynamics of the system are the result of a conservation of the properties of said system.

With the success of PINNs, more and more theoretical research on the integration of physical principles has been taking place. Namely, how do different formulations of a system’s equations can be used as prior knowledge to bias our model. The most common way in mechanics is to represent a system as its state-space. An alternative formulation can be done from the point of view of the energy of the system through the Hamiltonian formulation. For our state-space $ \mathbf{z} $ in the case of a MDOF formed by the pair $ \left(\mathbf{q},\mathbf{p}\right) $ , the Hamiltonian is formulated as:

(18) $$ H\left(\mathbf{q},\mathbf{p}\right)=\frac{1}{2}{\mathbf{p}}^T{\mathbf{M}}^{-1}\left(\mathbf{q}\right)\mathbf{p}+V\left(\mathbf{q}\right). $$

With $ {\mathbf{M}}^{-1}\left(\mathbf{q}\right)\succcurlyeq 0 $ . The Hamiltonian is considered to be separable when $ H\left(\mathbf{q},\mathbf{p}\right)=T\left(\mathbf{p}\right)+V\left(\mathbf{q}\right) $ . When looking at the variations of the Hamiltonian over time, one notices that $ \dot{\mathbf{q}}=\frac{\partial H}{\partial \mathbf{p}} $ and $ \dot{\mathbf{p}}=-\frac{\partial H}{\partial \mathbf{q}} $ . This leads to the Hamiltonian being time-invariant since:

(19) $$ \dot{H}={\left(\frac{\partial H}{\partial \mathbf{q}}\right)}^T\dot{\mathbf{q}}+{\left(\frac{\partial H}{\partial \mathbf{p}}\right)}^T\dot{\mathbf{p}}=0. $$

Greydanus et al. (Reference Greydanus, Dzamba and Yosinski2019) are the first to have introduced the Hamiltonian formalism to neural networks to bias the model for physical data. The method is closer to a Physics-informed than to a Physics-encoded model since the formalism is added through the loss function. The Hamiltonian loss is given by (we remind the reader that $ \Theta $ are the parameters of the model):

(20) $$ {\mathrm{\mathcal{L}}}_{\mathrm{HNN}}={\left\Vert \frac{\partial {H}_{\Theta}}{\partial \mathbf{p}}-\frac{\partial \mathbf{q}}{\partial t}\right\Vert}_2+{\left\Vert \frac{\partial {H}_{\Theta}}{\partial \mathbf{q}}-\frac{\partial \mathbf{p}}{\partial t}\right\Vert}_2. $$

In this case, is the neural network is learn the gradients of the system. Such a formalism prevents to the neural network prediction to stray away from the true state of the system by grounding it in the physical domain.

This property can also be found in symplectic integrators, that is, integrators derived from the Hamiltonian formalism. For a given initial value problem, a discrete integration of the system can be performed with explicit integrators such as Euler integration. Such integration relies on the local Taylor expansion of the flow of the system. For an integration of the nth order, an integration error of the n + 1-th order accumulates at every time step, leading to a drift from the true dynamics of the system. This drift can be avoid by opting for symplectic gradients instead. The first-order symplectic gradients can be derived as follows:

(21) $$ \dot{\mathbf{q}}=\frac{\partial H}{\partial \mathbf{p}}\Rightarrow {\dot{\mathbf{q}}}_t=\frac{{\mathbf{q}}_{k+1}-{\mathbf{q}}_k}{h}\Rightarrow {\mathbf{q}}_{k+1}={\mathbf{q}}_k+h\frac{\partial H}{\partial {\mathbf{p}}_k}. $$

(22) $$ \dot{\mathbf{p}}=-\frac{\partial H}{\partial \mathbf{q}}\Rightarrow {\dot{\mathbf{p}}}_t=\frac{{\mathbf{p}}_{k+1}-{\mathbf{p}}_k}{h}\Rightarrow {\mathbf{p}}_{k+1}={\mathbf{p}}_k-h\frac{\partial H}{\partial {\mathbf{q}}_k}. $$

For $ t $ the continuous time, $ k $ the discrete time, and $ h $ the discretisation. For the coefficients of symplectic integrators of the nth order (equivalent to their explicit counterparts), we refer to the formula derived by Yoshida (Reference Yoshida1990) ( $ {c}_i $ and $ {d}_i $ given in the paper):

(23) $$ {e}^{h\left(T\left(\mathbf{p}\right)+V\left(\mathbf{q}\right)\right)}=\prod \limits_{i=1}^n{e}^{c_i hT\left(\mathbf{p}\right)}{e}^{d_i hV\left(\mathbf{q}\right)}+\mathcal{O}\left({h}^{n+1}\right). $$

The symplectic integration principle has been extended to a plethora of physics-encoded models. Chen et al. (Reference Chen, Zhang, Arjovsky and Bottou2020) modified the original HNN by Greydanus by replacing the ANN with an RNN and by updating the gradients in a symplectic manner. Saemundsson et al. (Reference Saemundsson, Terenin, Hofmann and Deisenroth2020) train a Neural ODE with a split latent space with symplectic integrators. Sanchez-Gonzalez et al. (Reference Sanchez-Gonzalez, Bapst, Cranmer and Battaglia2019) estimate the Hamiltonian of the system, then derive the gradients from its estimate to update the model. David and Méhats (Reference David and Méhats2021) employ the same principle but replace the ANN with a Graph Neural Network (GNN). One of the issues with the Hamiltonian approach is that it assumes constant energy within the system, something that is usually not the case for most real-world applications. Many methods attempt to resolve this issue by incorporating dissipation into their formulation. Sosanya and Greydanus (Reference Sosanya and Greydanus2022), an adaptation of the original HNN that splits the gradients into their dissipative and non-dissipative components. Desai et al. (Reference Desai, Mattheakis, Sondak, Protopapas and Roberts2021) adopt the port-Hamiltonian formalism to adopt HNNs, making them apt to learn the dynamics of control systems.

The first work to extend this notion to DL is that of Greydanus et al. (Reference Greydanus, Dzamba and Yosinski2019), which enforces a symmetric gradient on a neural network trained to predict the dynamics of a conservative system. Saemundsson et al. (Reference Saemundsson, Terenin, Hofmann and Deisenroth2020) showed that the Hamiltonian formalism could be combined with the previously mentioned Neural ODEs, yielding the so-called Symplectic (state-space area preservation) Neural ODEs. Zhong et al. (Reference Zhong, Dey and Chakraborty2020) also used neural ODEs to learn the physics in an inferrable manner, applying Hamiltonian dynamics. Particularly, they parameterized the model in order to enforce Hamiltonian mechanics, even when only velocity data can be accessed as opposed to momentum. Bacsa et al. (Reference Bacsa, Lai, Liu, Todd and Chatzi2023) propose a method extending this reasoning to stochastic learning, where a symplectic encoder learns an energy-preserving latent representation of the system, and opens up new considerations for physics-embedded NN architectures.

We extend our SDOF oscillator example to other tasks to demonstrate the use of Neural ODEs. Neural ODEs are trainable forward models in that a neural network is used to approximate the flow of the system of study. In this context, the neural ODE estimates the latent space $ {\mathbf{z}}_t $ given the initial value problem starting at $ \left({\mathbf{z}}_0,t\right) $ . The Neural ODE flow estimation is done using the ResNet (He et al., Reference He, Zhang, Ren and Sun2015) such that the integration is accumulated on top of the residuals of the neural network. The neural network is optimized using the adjoint sensitivity method (Pontryagin et al., Reference Pontryagin, Mishchenko, Boltyanskii and Gamkrelidze1962). The results are seen in Figure 12.

Figure 12. Predicted versus exact solutions of state-space estimation of the Neural-ODE k + 1 predictor.

We can extend this method within the DMM framework, using the symplectic Neural ODEs, as per Bacsa et al. (Reference Bacsa, Lai, Liu, Todd and Chatzi2023). We change the problem accordingly: given that symplectic networks are made to deal with limit cycles, the forcing is switched from white excitation to a multisine excitation. The results are seen in Figure 13.

Figure 13. Predicted versus exact solutions of state-space estimation of the Symplectic Neural-ODE encoded DMM k + 1 predictor with uncertainty.