Hierarchical Validation Prevents Over-Fitting of the Neck Material Model for an Anthropomorphic Test Device Used in Underbody Blast Scenarios

In recent conflicts in Iraq and Afghanistan, 38% of all injuries were caused by improvised explosive devices (IEDs) [1]. IEDs create substantial vertical loads which can cause injury in the cervical spine and the head [1–4]. Protective equipment have increased survival rates, but these previously nonsurvivable injuries are now resulting in traumatic brain injury [4].

As such, having an anthropomorphic test device (ATD) that is biofidelic in underbody loading scenarios is important for evaluating injuries caused by underbody blasts. The most commonly used ATD is the hybrid III, which was developed for motor vehicle collisions [5], and as such, is generally not considered biofidelic in a vertical loading scenario [6–8].

In response to the need for an ATD that is biofidelic in an underbody loading scenario, the U.S. Army has created the Warrior Injury Assessment Manikin (WIAMan) program. The WIAMan is an ATD that is designed to be biofidelic in underbody loading conditions for the average soldier [9]. To ensure biofidelity, various cadaver tests have been performed to validate the ATD [10–13]. Three subinjurious under-body blast loading conditions were simulated with whole-body postmortem human surrogates (PMHS) and the WIAMan ATD. The PMHS were instrumented at 27 different locations and the response signals were evaluated against the ATD response using correlation analysis (CORA). CORA scores for the ATD neck ranged from 0.37 to 0.48 and scores for the head ranged from 0.54 to 0.75 indicating good biofidelity [10].

Various portions of the WIAMan have been characterized and implemented into the FEA model [14–16]. However, the head and neck (HN) has not yet been thoroughly characterized and validated against the ATD model. A properly modeled and validated neck material model is important because it influences metric values such as neck moment or force, as well as head accelerations. These metrics are used to ascertain injury risk and inaccurate models may result in inaccurate assessments of risk.

The ATD neck material is a butyl rubber with a measured hardness of 45 Shore A. An initial neck material model (Blatz-Ko rubber model, shear modulus G = 1.96 MPa and Poisson's ratio, ν = 0.463) was fit to data from testing on the molded neck part [17]. It was found to perform well in the constrained setup, but performed poorly in unconstrained setups. This combination indicates the model was over-trained in the neck-only (NO) data, similar to an over-fit mathematical model often described in the statistics field.

In this study, we aimed to improve the neck material model using a hierarchical approach to ensure accuracy and consistency across various scales. Material models were evaluated at the component level and then verified using untested, blinded component and whole body (WB) tests.

To improve the neck material model over the baseline Blatz-Ko rubber, we followed a hierarchical process: (1) created several potential material models based upon experimental material characterization data, (2) evaluated each potential material model at a component level and subsystem level using FEA, and (3) using the best performing material model from the component and subsystem level FEA, we validated its performance in several new conditions that were not used in the selection of the material model. This final step was intended as a blind validation of the selected material.

All FEA was performed using LS-Dyna R 8.0.0 (double precision). LS-Dyna is an explicit FEA program that is often used in highly transient events such as motor vehicle collision or underbody blasts.

Single element simulations were used to fit a material model to experimental data. Material fits were obtained using mcalibration software package (Veryst Engineering, Needham Heights, MA). mcalibration uses a variety of optimization techniques to minimize the normalized mean absolute difference between the single-element simulations and experimental data.

Material characterization on butyl rubber samples from a batch of the same rubber formulation was performed by Johns Hopkins University Applied Physics Laboratory (JHUAPL). The Bergstrom-Boyce model performed very well for this material but required slightly different parameters to achieve a proper fit for each mode of loading/unloading. Therefore, seven experimental tests were performed at strain rates ranging from 0.01/s to 100/s in tension and 0.01/s to 100/s in compression to determine the best parameters for each mode of loading and then a lookup table was generated. For the compression tests, the rubber specimen was cylindrical in shape (28.6 mm diameter × 12.5 length) and cyclic testing to 30% strain using a triangle wave for two cycles was used. The samples were placed between two lubricated platens and displacement was measured using the MTS stroke. Dog-bone samples were used for tension tests with 2 mm thickness × 6.25 mm width. The samples were centered in two clamps in glass filled nylon low mass pin joint fixtures and tested to failure. Local strain was determined from 1D DIC (digital image correlation) from high speed camera images and was then averaged over an area of interest. Engineering strain was measured on a local basis since the polymers experience large amounts of deformation.

The mesh of the WIAMan ATD neck model consists of 207,000 elements, and 195,000 nodes. All of the elements are solids and are comprised of 72% hex elements and 28% tetrahedral elements (around the face and skull).

Material fits were trained and validated on the component level (NO), and subsystem level (HN). The selected material was then validated on the subsystem level (two additional configurations) and the whole body level. With the exception of the whole body tests, each test was conducted at three different displacement or velocity speeds, labeled V1–V3. These configurations consisted of the following:

Neck-only setup: Fully constrained on both the top and bottom (Fig. 1(a)). The adaptor plates are modeled for placement between the test device neck and the load cells. The upper load cell was constrained, but the lower load cell was allowed upward (Z) translation. The input pulse is displacement and is applied to the bottom load cell. The displacement is controlled with a ramp (0.4, 0.6, and 0.8 mm/ms), followed by a hold (70 ms). Upper neck force and moment were recorded through cross section that were implemented in the ATD FEM model.

Head and neck setup: To mimic the experimental setup, the head and neck was placed on a tray system and a velocity is applied to the tray. The slider bearings on the linear bearing of the cart constrain the model to move only in the vertical direction. The neck was constrained on the bottom, but unconstrained on top (Fig. 1(b)). The input boundary condition is a velocity-controlled curve applied to the bottom of the plate as determined experimentally. The three speeds tested were 1, 2, and 3 m/s. For anatomical positioning, relative to the transverse plane, the nominal posture is flexed 3 deg.

Extended head and neck setup (HN-E): This setup is similar to the HN setup, except that the neck is extended 8.5 deg from the nominal position. This extension is measured about the neck adjustment plate in the model and is in agreement with previous experimental studies.

Flexed head and neck setup (HN-F): This setup is similar to the HN setup, except that the neck is flexed 19 deg from nominal position. This flexion is about the neck adjustment plate in the model and is based on previous experimental studies.

Whole body setup: The ATD is gravity settled onto a vertically accelerating seat. The seat is constrained to move only in the vertical (Z) mimicking underbody blasts and previous experimental setups. The ATD is fully unconstrained. This is most representative of actual underbody blast loading. The input condition is velocity controlled with separate pulses for the seat and the floor. Each of these pulses are determined from experimental conditions (Fig. 1(c)). Since the neck was the region of interest only responses in this area where examined.

Material selection was performed using the NO and HN setup (training datasets), and material model validation was performed using the HN-E, HN-F, and WB setups (validation datasets).

Simulation results were evaluated using CORA [18] without a corridor score and appropriate weights for signals according to their importance. All signal weights were decided a priori. The signal score is divided into three components: phase, size, and shape. Each of these components were weighted equally. The CORA score ranges from 0 to 1, with 1 being a perfect agreement between the simulation and experimental curves.

As CORA has shortcomings as an objective measure [19], primary peak characteristics were evaluated as well. Primary peak characteristics consisted of evaluating the peak magnitude and time-to-peak (TTP). For situations where there were multiple peaks, the primary peak was identified using the experimental data.

The NO setup was evaluated across six channels: lower neck force X, lower neck force Z, lower neck moment Y, upper neck force X, upper neck force Z, and upper neck moment Y. All the other setups were evaluated across seven channels: head velocity, head rotation Y, lower neck force Z, lower neck moment Y, lower neck rotation Y, upper neck force Z, and upper neck rotation Y. At a given speed, an average CORA, TTP, and Peak score were calculated across all signals. The combined total score is the sum of all speeds and was used for material selection. The selected material model was then evaluated in the validation of the model in the extended, flexed, and whole body condition. While results in the training data are presented against the baseline model, the validation results are only presented for the selected material model.

Five different Bergstrom Boyce material models were fit to various components of the experimental data (Tables 1 and 2). The Bergstrom Boyce is a nine parameter material model with viscoelastic and viscoplastic effects [20]. Because Bergstrom Boyce captured the predominant mechanical behaviors, notably the strain rate dependence and hysteresis, no other material models were fit. Individual results from the COL dataset from MCalibrate are shown below in Fig. 2. In general, there was good agreement between the experimental testing (solid lines) and the single element simulations (dashed lines).

Training simulation results (configurations NO and HN) are shown in Figs. 3 and 4 and are tabulated in the Appendix in Tables 4–6. Time history traces of the tests are found in the Appendix, in Figs. 5 and 7. Given the scoring system introduced, the maximum possible score for the training set was 9 (three different velocities, max score of 3 per velocity trial). Analyzing results from the neck-only condition (Fig. 3), the all-compressive (AC) material model had the best overall score with 7.70, closely followed by combined one loop (COL) material model with 7.53. The baseline Blatz-Ko material model performed second worst and had an overall score of 7.24, although, as noted in Fig. 3, the model preformed on par with the updated material fits due to the highly constrained boundary conditions.

Focusing on Fig. 4, however, which was a free condition of the head and neck, there is much clearer separation between the material model's performance over the baseline. Even though the all-compressive material model had the best overall score, the combined one loop and all tensile truncated (ATT) model preformed similarly. The combined one loop material model was selected for further validation because it captured both tension and compression responses. The all-compressive and all tensile truncated materials were only fit to compression response and tension response, respectively. Among the three highest scoring material fits, the difference in score was minimal. In this case the difference between the original Blatz-Ko (BK) (5.41) and combined one loop (6.72) material was more pronounced, this is also noted in the time history traces in the Appendix.

In the validation data (configurations HN-E, HN-F, and WB), the COL material model performed similarly to the training datasets. The CORA values ranged from 0.707 to 0.793, compared to 0.637 to 0.726 for the training dataset. The TTP scores ranged from 0.824 to 0.893, compared to 0.806 to 0.904 in the training dataset. The Peak scores ranged from 0.726 to 0.915 in the validation datasets, compared to 0.797 to 0.881 for the training datasets. These results are summarized in Table 3.

The baseline material selected for the model, a Blatz-Ko rubber formulation, was found to be limited in accurately predicting realistic physical response across multiple validation simulations. This deterioration of performance is evident between the NO and HN conditions with total scores of 7.25 versus 5.41 out of a possible total of 9. The Blatz-Ko formulation in LS-Dyna is relatively simple, with an assumed Poisson's ratio of 0.463, and is effectively a single parameter model based on material shear modulus. This is reasonable for the nearly incompressible rubber of the neck model; however, it was unable to capture the viscoelastic response of the material, as noted by the higher levels of oscillation observed in Fig. 6 where the condition was free. Incidentally, the experimental oscillation in Fig. 5 was thought to be related to specifics of the test setup that were not modeled. In our case, the model action was applied as a prescribed motion of the boundary plate based on experimental data.

As the C material model was fit to all experimental data, and, as such, contained the most range of possible stress–strain states, it was surprising that it did not perform well. It was also noteworthy that the material model AC, which was fit only to compression material tests, performed the best. We believed that this was due to a compressive loading bias in our training modes. This potential compressive loading bias led us to choose the COL material model, which included compressive and tensile material characterization data and performed similarly to the AC material model, to use as the WIAMan neck material model. Thus, indicating that while the combined score presented here is a robust means to evaluate candidate material models, engineering judgement can certainly play a role. In this case, the presence of the compressive data in the material fit outweighed the 0.05 difference in combined score (Table 6).

Generally, higher velocities saw better agreement with the experiment, and this was the most notable in the COL material model. This is reflected in the CORA scores (0.556–0.777), TTP scores (0.727–0.795), and peak scores (0.750 to 0.832). CORA had the largest increase; this was due to the oscillations observed in the simulations, but not in the experiment, that CORA penalized.

Good agreement with the experimental data reinforces our choice of a Bergstrom Boyce material model. The Bergstrom Boyce material model was chosen because it exhibits characteristics observed in the experimental data, notably viscoelastic and viscoplastic effects. Other material models generally do not inherently have these characteristics, or if they do, are often phenomenological models.

By utilizing the combined metric approach, we are able to take into account peak and time to peak into our combined score. CORA is a good assessment tool when looking at signals holistically (over the defined evaluation interval), but the software alone would not have specifically evaluated how well peaks or time to peaks match from experiment to simulation. Many injury metrics, especially those in the head, utilize peak accelerations when determining injury risk [21]. If CORA was used without the peak and TTP evaluation, there is the possibility of good overall curve agreement, without good peak or TTP agreement, which are more important when evaluating the risk of injury [22].

Similar performance was seen across the HN-E and HN-F setups. Both of these matched the HN setup, but with a change to neck angle. This may indicate that they are not substantially different and may not have been a good validation dataset. However, the WB still shows similar results, so we do not believe the model is over-trained.

The validation datasets were not used to evaluate the model beforehand, and were not considered during the material model selection. By doing so, we were able to evaluate the ability of the material model to extrapolate onto untrained scenarios. If we saw model performance drop significantly for the validation datasets, this would have indicated an over-trained material model.

Importantly, we measured similar scores between the training and validation datasets, which indicates that our model is not overly trained to a particular set of data. Additionally, like the training simulations, scores increased at higher velocities, and this further evidenced the score stability between the training and validation datasets. While we cannot say how the model will perform in very different scenarios than what are evaluated here, for instance, a side-impact loading condition, we expect scores for similar underbody loading scenarios would be similar to the scores the model received for the training and validation datasets. This is important because the model—and ATD itself –is intended to be placed into exploratory scenarios where it may not be feasible to validate. These scenarios will inform design decisions based on evaluated injury risk to the modern warfighter.

At low velocities, we saw oscillations that are indicative of an underdamped mass–spring–damper system. This is most likely due to other components within the head and neck that do not dissipate enough energy. This is evidenced by the good agreement in the neck-only setup, but once the head and neck are combined, we see these oscillations. Additionally, the oscillations were not observed in the whole body simulations, which are more representative of actual underbody loading scenarios. We envision that these oscillations can be addressed through various mechanisms. Furthermore, the lower neck rotation (L Neck Ry) is a limitation. While the highly transient initial behavior of the lower neck rotation is not captured (Fig. 6), the overall rotation in this case is small, less than half a degree, and the overall trend predicted by the FEA model is reasonable.

One mechanism is training the material model based upon more accurate strain rates experienced during testing. Material properties often vary significantly at different strain rates and tuning error can be reduced by training the material model to the strain rates expected during typical loading scenarios. Material characterization experimentation was performed at expected strain rates, but, in the absence of any superior material model, these expected strain rates were calculated based upon the underperforming initial Blatz-Ko material model. The five different Bergstrom–Boyce material models were selected from the previous experimental testing done, however, other Bergstrom–Boyce models were not explored and is a limitation of the study. However, by using a better material model, such as the COL, we can recalculate expected strain rates, and then redo material characterization experiments at those strain rates. Once the expected strain rates and the calculated strain rates from the best fit material match we can assume that our strain rate calculations are accurate.

This work focused on the validation of the head and neck kinematics of a finite element model of an anthropometric test device for use in underbody blast events. While the focus of the work was on the head and neck, the general practice is broadly applicable. By using a hierarchical approach with a combined scoring system, we achieved a robust model that performed similarly in validation datasets compared with training datasets. A combined scoring system that incorporates the well-known CORA evaluation system with explicit comparisons on signal peak and time to peak was introduced. The combined scoring system is unique and valuable because it captures both model response to experimental data (biofidelity) as well as signal peak (used in injury risk assessment). The model exhibited a 9% point increase in evaluation score versus the baseline Blatz-Ko material model. This leads to confidence that, within reason, the model can be extrapolated to new conditions. While we felt that the material model had shortcomings at lower velocities, it performed exceptionally well for higher velocities—which are the velocities typically encountered during an underbody blast scenario.

The head neck model geometry was developed by Humanetics ATD for the WIAMan Engineering Office of the U.S. Army Research Lab and consists of the head, neck, and personal protective equipment (PPE). The model of the experimental setup was developed by Duke University (neck-only). We would like to thank Alex Iwaskiw and Matt Shanaman (JHUAPL) for their work providing material characterization. We would like to thank Dr. Dale Bass and Maria Ortiz (Duke University) for performing the NO component tests. We would like to thank David Weyland and Hollie Pietsch (CCDC GVSC) for performing the head-neck ATD experimental tests (HN, HN-F, HN-E). We would also like to thank Dr. Dean Demetropoulos and Kyle Ott for performing the whole body experiments (WB). The authors gratefully acknowledge Jeremy Schap, formerly of Wake Forest University Center for Injury Biomechanics for his contribution to this work as well as the support of the administration of the Wake Forest University DEAC cluster, Adam Carlson and Damian Valles.

• US Army Research Lab via Johns Hopkins University Applied Physics Lab, Warrior Injury Assessment Manikin Development (Grant No. GTS: 40504; Funder ID: 10.13039/100012314).