The first part is the input data generation. For the conventional PINNs, points have to be generated 1) in the interior of domain to satisfy the PDE residual, 2) on the boundary of domain to satisfy the boundary conditions, and 3) additional points to satisfy the initial conditions. However, in the modified approach, since the output is adapted so as to satisfy all of the conditions simultaneously, only the interior points are required to be generated. The part of the code generating the interior data points by Latin hypercube sampling has been illustrated in the following snippet.

Next, the fully connected deep neural net architecture is constructed according to the user-defined number of layers “numLayers” and number of neurons per layer “numNeurons.” The trainable parameters (weights and biases) for every layer is initialized and stored in the fields of a structure array called “parameters.” The instance of initializing the weights and biases of the first fully connected layer has been captured by the following snippet. Here, the network weights are initialized by the He initialization [17] implemented by the function “initializeHe.” The He initializer samples the weights out of a normal distribution with zero mean and variance = 2 N i , where N i is the number of input channels. This function “initializeHe” takes in two input arguments, one is the size of trainable parameters “sz”and the other is “ N i ” and returns the sampled weights as a “dlarray” object. Note that “dlarray” is a built-in deep learning array in MATLAB employed for customizing the training process of DNNs. It enables numerous numerical operations including the computation of derivatives through AD. The network biases have been initialized by the zeros initialization implemented by the function “initializeZeros.” As it is evident that the initialization schemes can be easily customized, other initializers like Glorot or Xavier, Gaussian, orthogonal and others can be readily employed depending on the model type and choice of the user. In fact, a wide variety of initialization schemes of trainable parameters for various type of DNNs can be found in the MATLAB documentation. ¹

At this stage, the network is to be trained with user-specified value of parameters like, number of epochs, initial learning, decay rate along with several other tuning options. It is worth noting that multiple facilities to allocate hardware resources are available in MATLAB for training the network in an optimal computational cost. This include using CPU, GPU, multi GPU, parallel (local or remote) and cloud computing. The steps performed during the model training within the nested loops of epoch and iteration in mini-batches have been illustrated in the following snippet. To recall, an epoch is the full pass of the training algorithm over the entire training set and an iteration is one step of the gradient descent algorithm towards minimizing the loss function using a mini-batch. As it can be observed from the snippet that three operations are involved during the model training. These are 1) evaluating the model gradients and loss using “dlfeval” ² by calling the function “modelGradients” (which is explained in the next snippet), 2) updating the learning rate with every iteration and epoch and 3) finally updating the network parameters during the backpropagation using adaptive moment estimation (ADAM) [18]. In addition to ADAM, other stochastic gradient descent algorithms like, stochastic gradient descent with momentum (SGDM) and root mean square propagation (RMSProp) can be readily implemented via their built-in MATLAB routines.

The next snippet presents the function “modelGradients.” This sub-routine is the distinctive feature of PINNs where the physics of the problem is encoded in the loss functions. As mentioned previously, in conventional PINNs, the system response is assumed to be a DNN such that U=modelU(parameters,dlX,dlY,dlT). A difference to the expression of U can be observed in this snippet where the DNN output is modified based on the ICs and BCs. As obvious, this modification will change from problem to problem. In this case, the expression is shown for illustration and is related to Eq. 33 of Example 4 defined in the next section. As the name “dlgradient” suggests, it is used to compute the derivatives via AD. After evaluating the gradients, the loss term enforcing the PDE residual is computed.

As the modified DNN output ensures the satisfaction of ICs and BCs, only the loss term corresponding to PDE residual is necessary. Instead, if conventional PINNs was used, separate loss terms originating from the ICs and BCs would have to be added to the residual loss. Finally, the gradients of the combined loss w.r.t. the network parameters are computed and passed as the function output. These gradients are further used during backpropagation.

As obvious, there will be another loss term involved while solving an inverse problem which minimizes the discrepancy between the model prediction and the measured data. The parameter to be identified is updated as another additional hyperparameter of the DNN along with the network weights and biases. This can be easily implemented by adding the following line: c_update = parameters.(“fc” + numLayers).opt_param; and evaluating the PDE residual as f1 = c_update*(Uxx + Uyy) - Utt in Example 4. In doing so, note that c in line 22 of the snippet will be replaced by c_update.

The “modelU” function has been illustrated in the next snippet. Here, the fully connected deep neural network (FC-DNN) model is constructed as per the dimensionality of input and network parameters. In particular, the fully connect operations are performed via “fullyconnect.” This function uses the weighted sum to connect all the inputs to each output feature using the “weights,” and adds a “bias.” Sinusoidal activation function has been used here. The sub-routine returns the weighted output features as a dlarray “dlU” having the same underlying data type as the input “dlXYT.”

Once the PINNs model is trained, it can be used to predict on the test dataset. It is worth noting that the deep learning library of MATLAB is rich and consists of a diverse range of built-in functions, providing the users adequate choice and modelling freedom. In the next section, the performance of conventional and modified PINNs is accessed for solving four representative structural vibration problems, involving solution of ODE including multi-DOF systems, and PDE. In doing so, both forward and inverse problems have been addressed. Complete executable MATLAB codes of PINNs implementation for all the example problems can be found in the Supplementary Material.

The forced vibration of the spring-mass system can be expressed by u ̈ + ω n 2 u = f 0 ⁡ sin ⁡ ω t (12) where u, ü, ω_n, f_n, ω and t represent displacement, acceleration, natural frequency, forcing amplitude, forcing frequency and time, respectively. The initial conditions are u(t = 0) = 0 and ü (t = 0) = 0, where ü represents the velocity. The analytical solution to the above system is given by u t = f 0 ω n 2 − ω 2 sin ⁡ ω t − r ⁡ sin ω n t (13)

where, r = ω/ω_n is the frequency ratio.

As mentioned previously, in the realm of the PINNs framework, solution space (of the ODE, for this case) can be approximated by DNN such that u ̂ = N u ( t , θ ) , where the residual of ODE is evaluated with the help of AD. Essentially, this is an optimization problem which can be expressed as, arg min θ ∈ R d L θ ≔ ∂ 2 u ̂ ∂ t 2 + ω n 2 u ̂ − f 0 ⁡ sin ⁡ ω t 2 + u ̂ t = 0 2 + ∂ u ̂ ∂ t t = 0 2 (14) where, ‖ ◦ ‖ denotes ℓ 2 -norm. For the numerical illustration, it is assumed that ω = 3 , r = 1.5 and f 0 = 1 . The displacement u is approximated using a fully-connected neural network with 4 hidden layers and 20 neurons per layer. Sinusoidal activation function has been used due to the known periodic nature of the data [19]. 20,000 collocation points have been generated for time data t ∈ [ 0 4 π ] with the help of Latin hypercube sampling. The neural network is run for 1,000 epochs and the mini-batch size is 1,000. The initial learning rate is assumed to be 0.01 and the popular ADAM optimizer is employed. For testing the PINNs framework, 5,000 points were uniformly generated for time t ∈ [ 0 4 π ] . The solution u ̂ obtained using the PINNs framework has been compared with the actual (analytical) solution u in Figure 2A.

It can be observed from Figure 2A that the conventional PINNs framework is not capable of capturing the time response variation satisfactorily. As discussed in the previous sections, the reason is related to the regularization of the loss term in Eq. 14 and has been recently addressed in [16]. Although their approach proved to be effective, it entails extra computational effort.

Therefore, an alternative approach has been employed in this work to address the scaling issue which requires no additional computational cost compared to that of conventional PINNs. For avoiding multiple terms in the loss function, a simple scheme for modifying the neural network output has been adopted so that the initial and/or, boundary conditions are satisfied. To automatically satisfy the initial conditions in the above problem, the output of the neural network u ̂ is modified as, u ̂ mod = t u ̂ (15)

Since the modified neural network output is u ̂ mod , the new loss function can be expressed as, arg min θ ∈ R d L new θ ≔ ∂ 2 u ̂ mod ∂ t 2 + ω n 2 u ̂ mod − f 0 ⁡ sin ⁡ ω t 2 (16)

Following this approach, significant improvement in approximation of the displacement response has been achieved as shown in Figure 2B. Next, the implementation of PINNs has been illustrated for an inverse setting. For doing so, the same problem as defined by Eq. 12 is re-formulated such that the displacement time history is given in the form of measurements and the natural frequency ω n has to be identified. The optimization problem can be expressed as, arg min θ ∈ R d , ω ∈ R L θ , ω ≔ ∂ 2 u ̂ ∂ t 2 + ω n 2 u ̂ − f 0 ⁡ sin ⁡ ω t 2 + u ̂ − u * 2 (17) where, u * represents the measured displacement data. 15,000 collocation points have been generated for time data t ∈ [ 0 4 π ] with the help of Latin hypercube sampling. 2,500 displacement data points were used for artificially simulating the measurement data and 5 % uniform random noise was added. The architecture and the parameters of the neural network is the same as the previous case. The results have been presented in the form of convergence of the identified parameter ω in Figure 2C. The converged value of ω = 3.0 demonstrates exact match with the actual value. It is worth mentioning that the PINNs framework is inherently adapted to also provide the solution to the ODE along with the identified parameter in the inverse setup. This demonstrates that the PINNs framework can be easily adapted for solving forward and inverse problems in structural vibration.

The second example concerns a forced vibration of a damped spring-mass system and can be expressed by u ̈ + 2 ζ ω n u ̇ + ω n 2 u = f 0 ⁡ sin ⁡ ω t (18) where u, u ̇ , ü, ω_n , ζ, f₀, ω and t represent displacement, velocity, acceleration, natural frequency, damping ratio, forcing amplitude, forcing frequency and time, respectively. The initial conditions are u(t = 0) = 0 and u ̇ (t = 0) = 0. The analytical solution to the above system can be found in [20].

As mentioned previously, in the realm of the PINNs framework, solution space (of the ODE, for this case) can be approximated by DNN such that u ̂ = N u ( t , θ ) , where the residual of ODE is evaluated with the help of AD. Essentially, this is an optimization problem which can be expressed as, arg min θ ∈ R d L θ ≔ ∂ 2 u ̂ ∂ t 2 + 2 ζ ω n ∂ u ̂ ∂ t + ω n 2 u ̂ − f 0 ⁡ sin ⁡ ω t 2 + u ̂ t = 0 2 + ∂ u ̂ ∂ t t = 0 2 (19) where, ‖ ◦ ‖ denotes ℓ 2 -norm. For the numerical illustration, it is assumed that ω = 3 , ζ = 0.025 , r = 1.5 and f 0 = 1 . The displacement u is approximated using a fully-connected neural network with 4 hidden layers and 20 neurons per layer. Sinusoidal activation function has been used due to the known periodic nature of the data [19]. The neural network is run for 1,000 epochs and the mini-batch size is 1,000. The initial learning rate is assumed to be 0.01 and the popular ADAM optimizer is employed. The solution u ̂ obtained using the PINNs framework has been compared with the actual (analytical) solution u in Figure 3A. 20,000 collocation points have been generated for time data t ∈ [ 0 8 π ] with the help of Latin hypercube sampling to obtain the results in Figures 3A, B. For testing the PINNs framework, 5,000 points were uniformly generated for time t ∈ [ 0 8 π ] to obtain the results in Figures 3A, B.

It can be observed from Figure 3A that the conventional PINNs framework is not capable of capturing the time response variation satisfactorily. As discussed in the previous sections, the reason is related to the regularization of the loss term in Eq. 14. Therefore, to automatically satisfy the initial conditions, modified output of the neural network u ̂ mod is the same as Eq. 15 as the initial conditions are identical to that of the first example. Therefore, the new loss function can be expressed as, arg min θ ∈ R d L new θ ≔ ∂ 2 u ̂ mod ∂ t 2 + 2 ζ ω n ∂ u ̂ mod ∂ t + ω n 2 u ̂ mod − f 0 ⁡ sin ⁡ ω t 2 (20)

Following this approach, significant improvement in approximation of the displacement response has been achieved as shown in Figure 3B. The displacement response is presented over extended time in Figure 3C so as to investigate the performance of PINNs on the steady state response after the transients have died out. For generating the result in Figure 3C, 60,000 collocation points have been generated for the time data t ∈ [ 0 50 ] for training the network. For testing the PINNs framework, 40,000 points were uniformly generated for time t ∈ [ 0 50 ] . The approximation by PINNs is found to be excellent in terms of capturing the response trends.

Next, the implementation of PINNs has been illustrated for an inverse setting. For doing so, the same problem as defined by Eq. 18 is re-formulated such that the displacement time history is given in the form of measurements and both natural frequency ω n and damping ratio ζ have to be identified simultaneously. The optimization problem can be expressed as, arg min θ ∈ R d , ω ∈ R L θ , ω ≔ ∂ 2 u ̂ ∂ t 2 + 2 ζ ω n ∂ u ̂ ∂ t + ω n 2 u ̂ − f 0 ⁡ sin ⁡ ω t 2 + u ̂ − u * 2 (21) where, u * represents the measured displacement data. 10,000 collocation points have been generated for time data t ∈ [ 0 4 π ] with the help of Latin hypercube sampling. 1,000 displacement data points were used for artificially simulating the measurement data and 1 % uniform random noise was added. The architecture and the parameters of the neural network is the same as the previous case.

The results have been presented in the form of convergence of the identified parameters (natural frequency and damping ratio) in Figure 4. The converged value of ω n = 2.9985 and ζ = 0.0097 demonstrate close match with the actual values of 3 and 0.01, respectively. This demonstrates that the PINNs framework can be easily adapted for solving forward and inverse problems in structural vibration.

A 2-DOF lumped mass system as shown in Figure 5 is considered in this example [20]. This example has been included to illustrate the application of PINNs in a multi-output setting for the inference and identification of multi degree of freedom systems. The governing ODE and the initial conditions are as follows, m 1 0 0 m 2 q ̈ 1 q ̈ 2 + c 1 + c 2 − c 2 − c 2 c 2 + c 3 q ̇ 1 q ̇ 2 + k 1 + k 2 − k 2 − k 2 k 2 + k 3 q 1 q 2 = f 1 f 2 (22) with initial conditions q ( 0 ) and q ̇ ( 0 ) . In Eq. 22, q i , q ̇ i , q ̈ i represent the displacement, velocity, acceleration of the i th DOF, respectively, m i and f i represent the mass of the i th DOF and force acting at the i th DOF, respectively, c j and k j are the damping and stiffness coefficient of the j th connecting element, respectively. For this 2-DOF system, i = 1,2 and j = 1,2,3 . Since the free vibration problem has been undertaken, the right hand side of Eq. 22 is zero. Two cases of the free vibration problem have been considered, undamped and damped. For each of these cases, both forward and inverse formulations have been presented. The analytical solution to the above governing ODE considering undamped and damped cases, respectively, can be determined as, q t = ∑ i = 1 n d i u i ⁡ cos ω i t − ϕ i (23) q t = ∑ i = 1 n d i u i e − ζ i ω i t ⁡ cos ω d i t − ϕ i (24) where, constants d i and ϕ i have to be determined from the given initial conditions. n represents the number of DOFs, therefore n = 2 for the above system. ω i and u i are the i th undamped natural frequency and mode shape vector, respectively, obtained from the modal analysis. In Eq. 24, ζ i and ω d i represent the i th damping ratio and damped natural frequency, respectively.

As opposed to the previous examples, in general, the response associated with each DOF has to be represented by an output node of (multi-output) FC-DNN. Since the above example is a 2-DOF system, the response of the two DOFs are represented by two output nodes of an FC-DNN in the realm of PINNs architecture such that q ̂ = q ̂ 1 q ̂ 2 = N q ( t , θ ) . The optimization problem can be expressed as, arg min θ ∈ R d L θ ≔ M q ̂ ̈ + C q ̂ ̇ + K q ̂ 2 + q ̂ t = 0 − q 0 2 + q ̂ ̇ t = 0 2 (25) where, the gradients arising in Eq. 25 can be computed by AD. The following parameter values are adopted, m 1 = 9 , m 2 = 1 , k 1 = 24 , k 2 = 3 , k 3 = 0 , c 1 = 1 , c 2 = 0.125 , c 3 = 0 . An FC-DNN with 4 hidden layers and 20 neurons per layer is used. Sinusoidal activation function has been used due to the known periodic nature of the data. The neural network is run for 1,000 epochs and the mini-batch size is 1,000. The initial learning rate is assumed to be 0.01 and the popular ADAM optimizer is employed. Collocation points have been generated for time data t ∈ [ 0 8 π ] with the help of Latin hypercube sampling to obtain the results in Figures 6–8. For testing the conventional PINNs framework, 10,000 points were uniformly generated for time t ∈ [ 0 8 π ] to obtain the results in Figure 6. The undamped and damped time response obtained using conventional PINNs framework have been compared with the actual (analytical) solution in Figure 6.

It can be observed from Figure 6 that the conventional PINNs framework is capable of capturing the undamped and damped time response variation satisfactorily for two different ICs. The IC q 0 = q 01 q 02 = 1 3 , q ̇ 0 = 0 is adopted so that the system vibrates with the first natural frequency only as shown in Figures 6A, C, whereas the IC q 0 = q 01 q 02 = 1 0 , q ̇ 0 = 0 is a more general one resulting in a multi-frequency response as shown in Figures 6B, D. It is worth mentioning that the beat phenomenon exists in the free response of the above 2-DOF system due to close proximity of the two natural frequencies ( ω 1 = 2 and ω 2 = 2 ).

Next, PINNs has been implemented in an inverse setup for identification of system parameters both for the undamped and damped cases. For doing so, the same problem as defined by Eq. 22 is re-formulated such that the displacement time history data is available in the form of measurements and stiffness parameters ( k 1 and k 2 ) for the undamped case and stiffness and damping parameters ( k 1 , k 2 , c 1 and c 2 ) for the damped case, have to be identified simultaneously. The optimization problem for the undamped and damped case, respectively, can be expressed as, arg min θ ∈ R d , k 1 , k 2 ∈ R L θ , k 1 , k 2 ≔ M q ̂ ̈ + K q ̂ 2 + q ̂ − q * 2 (26) arg min θ ∈ R d , k 1 , k 2 , c 1 , c 2 ∈ R L θ , k 1 , k 2 , c 1 , c 2 ≔ M q ̂ ̈ + C q ̂ ̇ + K q ̂ 2 + q ̂ − q * 2 (27) where, q * represents the measured displacement data. Collocation points have been generated for time data t ∈ [ 0 8 π ] with the help of Latin hypercube sampling. 2,000 displacement data points were used for artificially simulating the measurement data and 1% uniform random noise was added. The architecture and the parameters of the neural network is the same as the forward formulation. The results have been presented in the form of convergence of identified system parameters in Figures 7, 8 for the undamped and damped case, respectively.

The converged values of k 1 = 24.0031 and k 2 = 2.9995 have been obtained from Figure 7 for the undamped case. The converged values of k 1 = 24.0064 , k 2 = 2.9995 , c 1 = 1.0030 and c 2 = 0.1248 have been obtained from Figure 8 for the damped case. The converged values of identified system parameters demonstrate close match with the actual values k 1 = 24 , k 2 = 3 , c 1 = 1 and c 2 = 0.125 . This demonstrates that the PINNs framework can be easily adapted for solving forward and inverse problems in multi-DOF systems. In addition to the adopted strategy to employ a single two-output FC-DNN to solve a 2-DOF system, two individual single output FC-DNNs were investigated. However, the latter failed to map the time response accurately due to the inability of two independent networks to adequately capture the dependencies of the coupled differential equations and hence, minimize the loss.

A rectangular membrane with unit dimensions excited by an initial displacement u = sin ⁡ π x ⁡ sin ⁡ π y has been considered in this example. The governing partial differential equation (PDE), initial and boundary conditions can be expressed as c ∇ 2 u = c ∂ 2 u ∂ x 2 + ∂ 2 u ∂ y 2 = ∂ 2 u ∂ t 2 ∀ x , y ∈ 0,1 , t > 0 (28) u = 0 ∀ x , y ∈ Γ u (29) u = sin ⁡ π x ⁡ sin ⁡ π y ∀ t = 0 , x , y ∈ 0,1 (30) ∂ u ∂ t = 0 ∀ t = 0 , x , y ∈ 0,1 (31) where, u is the displacement and c is the velocity of wave propagation. In Eqs 28–31, x , y represent the spatial coordinates, t represents time and Γ u denotes the spatial domain. The analytical solution to the governing PDE is u ( x , y , t ) = sin ⁡ π x ⁡ sin ⁡ π y ⁡ cos 2 π t .

Using the PINNs framework, solution of the PDE is approximated by a DNN such that u ̂ = N u ( x , y , t , θ ) , where residual of the PDE is evaluated with the help of AD. The optimization problem can be expressed as, arg min θ ∈ R d L θ ≔ c ∂ 2 u ̂ ∂ x 2 + ∂ 2 u ̂ ∂ y 2 − ∂ 2 u ̂ ∂ t 2 2 + u ̂ x , y ∈ Γ u 2 + u ̂ t = 0 − sin ⁡ π x ⁡ sin ⁡ π y 2 + ∂ u ̂ ∂ t t = 0 2 (32)

The displacement u is approximated using a fully-connected neural network with 4 hidden layers and 20 neurons per layer. Sinusoidal activation function has been used. 5,000 collocation points are generated for the spatial x , y ∈ [ 0 1 ] and temporal data t ∈ [ 0 1 / ( 2 2 ) ] with the help of Latin hypercube sampling. The neural network is run for 1,000 epochs and the mini-batch size is 1,000. The initial learning rate is assumed to be 0.01 and the popular ADAM optimizer is employed. For testing the PINNs framework, 1,000 points were uniformly generated for x , y and t . The solution in space obtained using the PINNs framework u ̂ (Figure 9B) has been compared with the actual (analytical) solution u (Figure 9A) for four different time instants t = 0.1 , 0.15 , 0.2 and 0.25.

It can be observed from Figure 9B that the conventional PINNs framework is not capable of capturing the time response variation satisfactorily. The reason is once again related to the regularization of the loss term in Eq. 32. The different terms related to the residual, initial and boundary conditions in the loss function are not satisfied simultaneously. Specifically, the fact that the condition u = 0 at the boundary of the domain not being satisfied in the predicted response by conventional PINNs can be visualized from Figure 9B.

To ensure the satisfaction of residual, initial and boundary conditions and improve upon the approximation accuracy, the neural network output has been modified as, u ̂ mod = t 2 x x − 1 y y − 1 u ̂ + sin ⁡ π x ⁡ sin ⁡ π y (33)

Since the modified neural network output is u ̂ mod , the new optimization problem can be expressed as, arg min θ ∈ R d L new θ ≔ c ∂ 2 u ̂ mod ∂ x 2 + ∂ 2 u ̂ mod ∂ y 2 − ∂ 2 u ̂ mod ∂ t 2 2 (34)

Following this modified PINNs approach, significant improvement in the spatial distribution of the displacement response has been achieved as shown in Figure 9C. Next, the implementation of PINNs has been illustrated in solving another inverse problem. For doing so, the same problem as defined by Eqs 28–31 is re-formulated such that the displacement time history is given in the form of measurements and the wave velocity c has to be identified. The optimization problem can be expressed as, arg min θ ∈ R d , c ∈ R L θ , c ≔ c ∂ 2 u ̂ ∂ x 2 + ∂ 2 u ̂ ∂ y 2 − ∂ 2 u ̂ ∂ t 2 2 + u ̂ − u * 2 (35) where, u * represents the measured displacement data. 25,000 collocation points have been generated for spatial coordinates x , y ∈ [ 0 1 ] and time t ∈ [ 0 1 / ( 2 2 ) ] with the help of Latin hypercube sampling. 5,000 displacement data points were used for artificially simulating the measurement data and 2 % uniform random noise was added. The architecture and the parameters of the neural network is the same as for the forward problem. The results have been presented in the form of convergence of the identified parameter c at the end of 10,000 epochs in Figure 9D. The converged value of c = 0.9902 demonstrates good match with the actual value c = 1.0 . It is worth noting that the PINNs framework is inherently adapted to also provide the solution to the PDE along with the identified parameter in the inverse setup. This demonstrates that the PINNs framework can be easily adapted for solving forward and inverse problems in structural vibration.