Fifty Shades of Brain: A Review on the Mechanical Testing and Modeling of Brain Tissue

Brain tissue is not only one of the most important but also the most complex and compliant tissue in the human body. While long underestimated, increasing evidence confirms that mechanics plays a critical role in modulating brain function and dysfunction. Computational simulations–based on the field equations of nonlinear continuum mechanics–can provide important insights into the underlying mechanisms of brain injury and disease that go beyond the possibilities of traditional diagnostic tools. Realistic numerical predictions, however, require mechanical models that are capable of capturing the complex and unique characteristics of this ultrasoft, heterogeneous, and active tissue. In recent years, contradictory experimental results have caused confusion and hindered rapid progress. In this review, we carefully assess the challenges associated with brain tissue testing and modeling, and work out the most important characteristics of brain tissue behavior on different length and time scales. Depending on the application of interest, we propose appropriate mechanical modeling approaches that are as complex as necessary but as simple as possible. This comprehensive review will, on the one hand, stimulate the design of new experiments and, on the other hand, guide the selection of appropriate constitutive models for specific applications. Mechanical models that capture the complex behavior of nervous tissues and are accurately calibrated with reliable and comprehensive experimental data are key to performing reliable predictive simulations. Ultimately, mathematical modeling and computational simulations of the brain are useful for both biomedical and clinical communities, and cover a wide range of applications ranging from predicting disease progression and estimating injury risk to planning surgical procedures.


Introduction
Brain tissue is one of the most complex tissues in the human body. Neurological disorders, including stroke, encephalitis, dementias, and epilepsy, have been identified as one of the major public health concerns by the world health organization. An additional threat are the consequences of neurotrauma with over 2 million people affected by traumatic brain injury each year [111]. While over the past several decades, neuroscience has mostly been limited to electrophysiological, biochemical, molecular and genetic studies [162], more recent studies show that mechanics plays a critical role for neuronal function and dysfunction [9,68]. External mechanical loads may cause damage of brain tissue during traumatic brain or spinal cord injury [101], but mechanical cues may also control developmental processes and progression of disease through mechanosensation of nerve cells [96,162]. For instance, mechanics seem to not only drive cortical folding during mammalian brain development [21,26,141,143], but also tumor cell migration and cell apoptosis [162]. Unphysiological changes in the 1 3 mechanical environment induce cortical malformations [21,27] associated with schizophrenia [76], autism [129], and epilepsy [14], or cell death in Alzheimer's patients associated with dementia [15].
Computational modeling and personalized simulations can provide fundamental insights into the underlying mechanisms during injury and disease. Such predictive models not only reduce the necessity of experiments on humans and animals, but also allow for the development of novel treatment strategies or the detailed planning of surgical procedures [171]. However, realistic predictions of mechanobiological processes in the brain require sophisticated mechanical models that capture the complex and unique characteristics of this ultrasoft, highly adaptive and heterogeneous tissue. While great efforts have been made to mechanically model the behavior of brain tissue in health and disease [68], contradictory experimental results have constantly hindered progress [29] and caused confusion and delay [7]. This review aims to explain these discrepancies and provide a comprehensive overview on the mechanics of brain tissue.
After addressing the challenges associated with brain tissue testing and modeling in Sect. 2, we review experimental observations based on different testing techniques to work out the diverse characteristics of brain tissue behavior under different loading conditions in Sect. 3. Subsequently, we present mechanical modeling approaches that are capable of mathematically describing the characteristic tissue behavior in Sect. 4. Finally, we discuss application-specific features to propose appropriate modeling approaches together with the corresponding set of material parameters in Sect. 5. We will close this review with a short conclusion and open challenges in Sect. 6.

Brain Tissue is Ultrasoft
With a shear modulus on the order of one kilopascal [22], brain tissue is ultrasoft-softer than any other tissue in the human body. For comparison, Fig. 1 highlights the stiffnesses of various organs in our body. Traditional mechanical test setups have originally been designed for much stiffer materials, even stiffer than bone. The extreme compliance of brain tissue pushes mechanical testing and modeling approaches to their limits. On the one hand, boundary conditions might be inappropriate. On the other hand, force sensors reach their sensitivity limit. As a result, early studies on the stiffness of brain tissue yielded significantly higher values than more recent studies [59]-the testing devices used at that time were merely not capable of recording the accurate tissue responses reliably.
To complicate matters, brain tissue deforms noticeably due to gravity-even just under its own weight. This could have important consequences on testing results, especially when the height of the specimen is relatively large [99]. It is thus difficult to control specimen geometry and local deformation states during biomechanical testing, which, in turn, will affect the forces recorded during testing [138]. This implies that the common assumption of homogeneous deformation during uniaxial tension or compression  [22,121,137,140] is debatable [147]. The only alternative to circumvent this concern are inverse studies, where finite element simulations are performed to capture local inhomogeneities in the deformation. However, it will always be difficult to exactly match the specimen geometry of experiment and simulation throughout the entire loading history. In general, ultrasoft tissues undergo large deformations-even under their own weight-and it is essential to formulate mechanical models with the framework of nonlinear continuum mechanics.

Brain Tissue is Highly Fragile
Due to the fragile nature of brain tissue, it is difficult to control tissue damage during biomechanical testing, especially when applying large deformations. In white matter brain tissue, axons, the long and slender nerve fibers, will tear when loaded above a certain threshold, as illustrated in Fig. 2 of [58]. Previously, this threshold has been determined at approximately 18% tensile strain. Measurements above, or even close to, this damage threshold are likely to render meaningless results.
In addition, the effect of adhesion, as illustrated in Fig. 3, enhances the risk of damaging specimens during preparation and handling [61]. This demonstrates that biomechanical tests need to be well-designed and performed with exceptional caution.

Brain Tissue is Biphasic
Unlike other soft tissues, the brain is not only ultrasoft, but also has an exceptionally high water content, 0.83 g/ml in gray matter and 0.71 g/ml in white matter [180]. From a total of about 80% water, approximately 20-40% is free flowing cerebrospinal fluid, while the rest resides inside the cells. Due to the biphasic nature of brain tissue, it is difficult to control drainage conditions during biomechanical testing [58]. As a result, depending on the testing setup, the incompressible fluid is trapped within the tissue or free to escape, which will largely affect the recorded response. This interaction between solid and fluid phases is highly sensitive to the time and length scales [29]. Depending on the type of loading and the testing setup, e.g., unconfined compression/extension versus nanoindentation experiments, different tissue components control the recorded mechanical stiffness. Nanoindentation experiments on intact brain tissue slices [20] closely mimic the nature of confined compression tests, where the fluid phase is mostly kept within the solid matrix and contributes to the tissue stiffness. In contrast, during unconfined compression and tension tests, the fluid is free to escape, and we primarily probe the elastic properties of the solid skeleton consisting of cells, intercellular connections, and extracellular matrix. Naturally, these tests introduce an additional time scale, the time for the fluid to escape the solid. Therefore, the results are sensitive to the type of preconditioning and the recovery time between two consecutive tests [23].

Brain Tissue is Highly Heterogeneous
Brain tissue displays a high degree of microstructural heterogeneity, as shown in Fig. 4. Its cellular and extracellular structure conforms with regionally varying functional demands, and so do its mechanical properties [173]. For most practical purposes, it seems too simplistic to consider brain tissue as a single homogeneous material. Figure 5 illustrates the major cell types in white and gray matter tissue. Gray matter regions contain mainly neuronal cell bodies, protoplasmic astrocytes providing neurons with nutrients, and microglia as active immune defense. White matter regions, on the contrary, contain axons, oligodendrocytes which wrap isolating myelin sheath around the axons, fibrous astrocytes, and microglia [27]. Notably, not only the cellular composition may be relevant for macroscopic tissue mechanics but also extracellular matrix components. The latter may show regional trends that even differ from those of brain cells.
Notably, the microstructure of brain tissue not only varies in space but also in time. When we learn a new task, for example juggling, the neurons responsible for fine motor skills will form further connections with one another and Fig. 2 Damage evolution and fracture process in a (14 mm/9.5 mm initial height/edge) prismatic specimen of white matter harvested from the occipital lobe in the frontal direction. Adapted from [58] Fig. 3 Illustration of adhesion effects during handling of brain tissue specimens create new synapses. Brain plasticity, the temporal changes in the tissue's microstructure, also affects the mechanical response. This emphasizes the complexity of brain tissue behavior and the high demands when aiming to develop mechanical models for brain tissue.

Brain Tissue is Not Easily Available
One major issue considering the mechanical testing of brain tissue is that it is exceptionally challenging to obtain adequate material. While it would be desirable to test the mechanical properties of brain tissue in vivoin its natural environment-there are major issues with in vivo testing [97]. Firstly, traditional mechanical testing setups, initially developed and used for non-organic materials, require to prepare specimens of a certain geometry-optimized to allow for homogeneous deformation states during testing [63]. This is obviously not possible without removing the material from its natural environment and destroying its biological functionality. Secondly, the calibration of mechanical models that are useful to study injury and disease require a versatile set Fig. 4 Microstructure of different brain regions, the cortex (C), the basal ganglia (BG), the corona radiata (CR), and the corpus callosum (CC). Klüver-Barrera (K.B.) staining was used to color myelin wrapped around nerve fibers in blue, neuropil in pink, and nerve cells in purple. Images are shown at 20 magnification. Adapted from [22] Fig. 5 Cellular organization of brain tissue. In white matter (corona radiata), myelinated axons allow for rapid nerve impulse conduction; intermediate oligodendrocytes connect and form several myelin sheaths. Fibrous astrocytes ensure supply of nutrients and synaptic processing. In gray matter (cortex), neurons form synapses with each other and with protoplasmic astrocytes. In both white and gray matter, microglial cells contribute to clearance of debris and synapse remodeling. Reprinted from [27] under Creative Commons Attribution License (CC BY) of different loading conditions, including deformations that go beyond physiological conditions and might induce tissue damage. Again, those tests can not be performed in vivo on a healthy individual.
On the other hand, the deformations imposed by common ex vivo testing techniques such as nanoindentation, uniaxial compression and extension, or simple shear will undoubtedly differ from those experienced in vivo. Recently, a newly established testing technique, magnetic resonance elastography, has enabled to study the mechanical properties of human brain tissue in vivo [112]. While magnetic resonance elastography is a powerful tool to perform non-invasive longitudinal studies to monitor the stiffness evolution under diseased conditions, for example in multiple sclerosis [182], it can not yet be immediately used to derive mechanical models for brain tissue, and fully replace ex vivo experiments [12]. Magnetic resonance elastography is based on the assumption of linear elasticity or viscoelasticity [146], and does not offer enough information to establish nonlinear, finite strain constitutive models for realistic computational simulations. Besides, the mechanical properties obtained from magnetic resonance elastography are sensitive to positioning, to the method to extract elastic parameters, and to the excitation frequency [177]. It is not at all trivial to accurately quantify regional mechanical properties, especially considering thin structures such as the cortex and the corpus callosum [131]. In summary, the better we understand the different "shades" of brain tissue behavior, the easier it will be to establish new appropriate testing techniques in the future.

Experimental Observations
From all these challenges, it becomes clear why experimental studies on the mechanical properties of brain tissue have been equivocal and why stiffness values in the literature have varied by an order of magnitude or more [119]. Depending on the testing techniques different aspects of the complex mechanics of brain tissue control the recorded response as illustrated in Fig. 6. In this section, we summarize the experimentally observed characteristics of brain tissue, discuss how different mechanical testing techniques may have resulted in contradictory findings, and critically evaluate facts that are still under debate. Ultimately, this will allow us to work out the key characteristics to establish appropriate constitutive models for brain tissue in Sect. 4. Figure 7 illustrates the loading-mode-specific elastic behavior of brain tissue-the response we expect for very slow processes such as tissue growth or progress of neurodegenerative diseases-here tested during unconfined compression, tension, and simple shear experiments for four different brain regions, the cortex, the basal ganglia, the corona radiata, and the corpus callosum as introduced in Fig. 4 [22]. Independent of the loading mode, we observe a strainstiffening behavior-the stiffness increases with increasing strain-which highlights the strong nonlinearity of the tissue response.   Figure 7 reveals that stresses are significantly higher in compression than in tension. This important characteristic of brain tissue has been consistently observed by several research groups [22,58,85,120,121]. A convenient test to capture the loading-mode-specificity of the mechanical response of brain tissue is to perform tests under combined loading conditions, as illustrated in Fig. 8 [22,133]: For combined compression/tension-shear experiments, the shear stresses increase noticeably with increasing compression, but only marginally with increasing tension [22]. We may attribute this behavior to the biphasic nature of brain tissue discussed in Sect. 2.3: Cerebrospinal fluid trapped inside the solid network of cells and extracellular matrix offers noticeably resistance in compression but only marginal resistance to tension. Figure 9 demonstrates that even for a relatively slow loading speed of 2mm/min corresponding to strain rates of approximately 0.0067 1/s, brain tissue exhibits a highly hysteretic response [22]. We expect this hysteresis to vanish for much slower rates, which might occur in vivo during brain growth or disease, but are difficult to reproduce during ex vivo testing. Figure 10 shows stress relaxation experiments for all four brain regions [22], which confirm the extreme timedependence of brain tissue with a stress relaxation of up to 80% within only 300 s. Interestingly, white matter tissue, with a stress relaxation of more than 70% after 300 s, is more viscous and responds more slowly than gray matter, with a stress relaxation of approximately 65%. Within white matter, specimens from the corpus callosum relax faster than specimens from the corona radiata. Within gray matter, the basal ganglia and the cortex exhibit a similar relaxation behavior. Stress relaxation percentages are slightly higher in shear than in compression, but both loading modes show similar regional dependencies. Average elastic response of brain tissue with standard deviations during multiple loading modes, simple shear, compression, and tension, for different brain regions: cortex, basal ganglia, corona radiata, and corpus callosum. Adapted from [22] 0.  9 Average unconditioned response of brain tissue with standard deviations during multiple loading modes, simple shear, compression, and tension, for different brain regions: cortex, basal ganglia, corona radiata, and corpus callosum. Adapted from [22] Fig. 10 Stress relaxation behavior of brain tissue with standard deviations and stress relaxation percentages after 5 minutes for simple shear and compression loadings in different brain regions: cortex, basal ganglia, corona radiata, and corpus callosum. Adapted from [22] 1 3

Brain Stiffness Increases with Increasing Strain Rate
Figures 11 and 12 illustrate that brain tissue not only stiffens with increasing strain but also with increasing strain rate. Rate-dependence of brain stiffness has consistently been reported in the literature using different testing techniques: shear testing or oscillatory shear testing [44,127,151], uniaxial compression or tension [57,83,85,132,137,140], indentation [20,107,136,163], and magnetic resonance elastography [36,145]. Figure 11 illustrates uniaxial compression, tension, and simple shear experiments performed at strain rates of 0.33 and 0.0067 1/s, respectively. Figure 12 shows nanoindentation experiments over a loading rate spectrum from 1 to 160 m/s . Within the analyzed loading rate regime both gray and white matter double their maximum forces and corresponding moduli when increasing the loading rate by two orders of magnitude [20]. This effect becomes particularly important for applications such as blunt or traumatic brain injury, where even higher strain rates occur [137,158]. Our results in Fig. 11 suggest that the effect of strain rate is more pronounced in white matter regions than in gray matter regions. This observation agrees well with the stress Fig. 11 Average initial loading segment with standard deviations for simple shear as well as uniaxial compression experiments loaded at different speeds of v = 100 mm/min and v = 2 mm/min correspond-ing to strain rates of approximately 0.33 and 0.0067 1/s, respectively, in four different brain regions: cortex, basal ganglia, corona radiata, and corpus callosum

Fig. 12
Loading rate sensitivity of gray and white matter. Sensitivity of indentation force versus indentation depth for varying loading rates reveals the rate-dependent nature of brain tissue. Indentation force and modulus increase with increasing loading rate. Adapted from [20] relaxation experiments in Fig. 10, which show that white matter relaxes faster than gray matter. A possible explanation for this behavior could be the difference in the permeability of gray and white matter. While for slower loadings, fluid has more time to escape, for faster loadings, it offers resistance, which leads to higher stresses. As white matter specimens seem to loose a higher amount of fluid during unconfined experiments than gray matter specimens [22,23,25], their strain rate effect is more pronounced. Figure 13 shows the conditioned tissue response during the third cycle of cyclic loading in simple shear, uniaxial compression, and uniaxial tension for four different brain regions, the cortex, the basal ganglia, the corona radiata, and the corpus callosum. Notably, during the third loading cycle, the stresses are significantly lower than during the first loading cycle, the unconditioned tissue response in Fig. 9. We conclude that brain tissue softens upon preconditioning. This observation agrees well with reported results, where the equilibrium shear modulus during the fifth cycle was consistently about 20% lower than during the first cycle [62]. Figure 14 illustrates that a characteristic softening between the first loading cycle and all subsequent cycles occurs whenever the tissue is loaded past the previous maximum strain. We once more attribute this characteristic conditioning to fluid leaking out of the sample: As the maximum strain increases, more fluid escapes the sample and no longer offers resistance to loading. In turn, the offset at zero strain upon unloading increases.

Brain Tissue Softens upon Preconditioning
For simple shear and tension, we observe jagged peaks during the first loading cycle of each strain level. We suggest two possible explanations for this phenomenon according to the challenges introduced in Sects. 2.1-2.4. Firstly, it is possible that part of the water molecules of the cerebrospinal fluid are non-covalently bonded to the highly charged macromolecules in the extracellular matrix; during initial loading, those weak interactions break and the water is free to escape. Secondly, these peaks could result from tissue that initially adheres to the upper or lower specimen holder, without being glued, but comes loose during loading.

Fig. 13
Average conditioned response of brain tissue with standard deviations during multiple loading modes, simple shear, compression, and tension, for different brain regions: cortex, basal ganglia, corona radiata, and corpus callosum. Adapted from [22] Figure 15 demonstrates that the characteristic conditioning behavior discussed in Sect. 3.5 is fully recoverable. When we load the same specimen after a recovery period of one hour, the cyclic response of brain tissue up to 10% compression follows almost the identical path of the initial loading-again showing the evident softening after the first loading cycle. A similar behavior has been observed during nanoindentation experiments: The tissue continuously softens throughout multiple indents, but fully recovers within a few hours of rest [20]. We conclude that the softening due to conditioning is not associated with tissue damage but with visco-or poro-elastic effects that recede over time. This is an important finding towards developing appropriate mechanical models for brain tissue behavior in Sect. 4. Notably, the tissue only recovers when it is kept hydrated. If the tissue was left in air, it would dry out and the response would stiffen over time [44].

Brain Tissue is not Notably Anisotropic
From a mechanics point of view, a fundamental question that needs to be answered is whether brain tissue is an isotropic Fig. 14 Representative stress versus amount of shear/stretch behavior in simple shear, compression, and tension for stepwise increasing shear/stretch-three cycles per shear/stretch level. The tissue exhibits substantial pre-conditioning during the first cycle of each level, and successively softens with increasing maximum shear/stretch or anisotropic material. Due to controversies regarding previous results in the literature, where some suggested that brain tissue was isotropic [22,22,91,120,137], but others showed that there were significant directional trends [52,53,85,134,165], we combined biomechanical testing of the same specimen in three orthogonal loading directions with antecedent diffusion tensor imaging to carefully analyze to which extent the fibrous microstructure of axonal networks in white matter results in an anisotropic macroscopic mechanical response [22]. This allows us to minimize falsifying effects of non-uniform fiber distributions or inter-specimen variations that might have affected previous studies. We determined nerve fiber distributions prior to specimen extraction and tested the same microstructurally anisotropic sample in multiple loading directions-the three possible modes for simple shear (FT, TT, and TF) and the two possible modes for compression and tension (FF and TT), as denoted in Fig. 16a.
While diffusion tensor images showed highly anisotropic diffusion properties, mechanical testing revealed no significant directional dependencies, as illustrated in Fig. 16. According to these results, while structurally anisotropic, brain tissue can be considered as an isotropic material from a mechanical point of view. This finding seems reasonable considering the fact that nerve fibers are primarily functional elements of our brain that only marginally contribute to its mechanical strength. We note, however, that mean stresses along fibers appeared slightly higher than perpendicular to the fibers under tensile loading and slightly lower under compressive loading. This could imply that-although statistically insignificant-nerve fibers contribute to tissue strength in tension but not in compression.
Our results agree well with an early study on the rheological shear response of human brain tissue [151], which showed higher directional variation in gray matter than in white matter, but neither of the differences appeared significant. Contradictory to our findings, studies on porcine brain tissue found a significantly stiffer shear response orthogonal to nerve fibers than along fibers in the corpus callosum [134]. In the corona radiata, however, the trend was opposite. The authors of this study sheared  where I denotes the normal vector of the face that is shifted by simple shear and J denotes the direction of shear, along (F) or transverse (T) to the fiber direction f 0 . c Direction-dependent behavior of brain tissue. Column plots of shear moduli (mean ± SD) averaged over three possible modes in simple shear (FT, TT, and TF), and two possible modes in compression and tension (FF and TT) for all specimens from the corpus callosum. The data indicate a marginally softer response along nerve fibers in compression and marginally stiffer response along fibers in tension but the differences appeared to be statistically not significant (n.s.). Reprinted with permission from [22] each specimen in two orthogonal directions similar to our experiments, but only estimated fiber orientations from anatomical knowledge and used rectangular specimen dimensions of 10 × 5 × 1 mm 3 . Notably, the measured shear stresses were consistently higher in the direction of the longer axis corresponding to the direction orthogonal to fibers in the corpus callosum, and to the fiber direction in the corona radiata. This could indicate that directional dependencies are an artifact of the non-cuboidal specimen dimensions rather than a result of the anisotropic distribution of nerve fibers, which could explain the contradictory results. Interestingly, yet other studies on the porcine corpus callosum found opposite trends with a significantly stiffer response in the fiber direction than perpendicular to it in dynamic shear [53] and tensile tests [165]. In both studies, specimens were relatively large with dimensions of up to 16 × 16 × 3 mm 3 and 5 × 5 × 60 mm 3 . Our diffusion tensor images showed that even in the much larger human brain, it would be challenging to extract specimens of that size that exhibit a sufficiently uniform microstructure. We therefore interpret the corresponding results with a degree of caution.
Overall, compared to other soft biological tissues including muscle, tendon, or arteries with pronounced directional stiffness differences of an order of magnitude or more, the directional stiffness of brain tissue only varies by a few percent and can reasonably well be assumed as isotropic. Nonetheless, while isotropy seems to be a valid assumption for the elasticity of brain tissue, diffusion or permeability properties might still be anisotropic. Furthermore, the role of axons could be decisive for white matter brain damage, where rupture of axons can lead to loss of brain function [37,66].

Brain Tissue Stiffness is Region-Dependent
Due to its high microstructural heterogeneity, brain tissue can hardly be considered as a single material with unified material properties. While early studies on the mechanical properties of the brain focused on brain tissue as a whole, more recent experimental studies have distinguished between different regions, i.e., white and gray matter [20,57,114,163], or, even more refined, cortex, basal ganglia, corona radiata and corpus callosum, as depicted in Fig. 4 [22,85]. Others tested the individual properties of the cerebrum, cerebellum, pons, and medulla [109].
Astonishingly, experimental results towards regional trends in brain tissue properties have been contradictory. While indentation experiments including those shown in Fig. 17 yield higher stiffnesses in white matter from the corona radiata than in the cortex [20,92,163], simple shear, unconfined compression and tension tests in Figs. 7, 8, 9, 10 and 11 show the opposite trend [22]. The origin of these, at first sight, contradictory results lies in the complexity of brain tissue behavior, which highly depends on the length and time scales of the applied loading. Furthermore, not only the tissue's elasticity shows local variations but also other properties such as viscosity and porosity. To disclose regional differences in brain tissue properties, it is therefore indispensable to thoroughly understand the different "shades" of brain tissue mechanics and their loading-mode specific effects.

Regional Trends Depend on the Loading Rate
Firstly, the rheological difference discussed in Sect. 3.4 leads to a rate-dependency of regional trends: White matter stiffens relative to gray matter with increasing loading rate. As a Fig. 17 Regional variation of gray and white matter moduli. Measurements at three different slices and three different locations reveal that the specimen moduli vary markedly across the brain. Gray matter, left, is softer than white matter, right, and displays smaller regional variations. Black horizontal lines indicate the mean; gray zones indicate the standard deviation. Adapted from [20] result, shear, compression or tension experiments in the fast loading regime associated with phenomena on the order of seconds or milliseconds such as traumatic brain injury have reported a higher stiffness for white matter from the corona radiata than for cortical gray matter [85,114], while experiments in the slow loading regime, as those shown in Fig. 11, indicate the reversed relationship. The strain-rate dependence of regional trends is also supported by a recent study showing that inter-regional mechanical properties become increasingly heterogeneous with increasing strain rate [108].

Regional Trends Depend on Drainage Conditions
To explain differences between nanoindentation and unconfined macroscopic experiments-both performed at relatively low rates-we refer to the biphasic nature of brain tissue discussed in Sect. 2.3. Nanoindentation of intact tissue slices can be regarded as a virtually confined test, where the free flowing fluid is mostly kept within the solid matrix. During simple shear, unconfined compression and tension tests, in contrast, the fluid is free to escape. As a result, unconfined tests reveal a markedly stiffer response for the unconditioned response, see Fig. 9, than for conditioned response, see Fig. 13, which emphasizes that the fluid plays a significant role in the tissue's mechanical strength: Brain tissue softens when part of the fluid escapes. In turn, when tissue samples have time to saturate themselves again after loading during a recovery period of one hour, in Fig. 15, they regain their initial stiffness.
When comparing the porous nature of different brain regions, our experiments indicate that the largest amount of fluid escapes from white matter specimens of the corpus callosum, closely followed by the corona radiata, while gray matter specimens from the cortex loose the least amount of fluid [23]. This difference in the permeability of gray and white matter is also reflected in larger hysteresis areas in white matter than in gray matter during cyclic loading [22]. We can rationalize these observations with the underlying tissue microstructure: the corpus callosum consists of a sparsely cross-linked network of unidirectional fibers, whereas the cortex consists of a densely connected network of dendrites that traps the fluid phase inside the tissue.
With regard to these observations, we suppose that during unconfined experiments, we primarily probe the "elastic" properties of the solid skeleton, the cells and the extracellular matrix. Those are stiffer in gray matter from the cortex than in white matter from the corona radiata, see Fig. 9 [22]. We observe a different behavior for confined compression experiments during nanoindentation with relatively large indenters, as shown in Fig. 18: The fluid cannot escape freely and we probe both solid and fluid phases. Due to the differences in the porous nature of gray and white matter, we record larger stiffnesses for white matter from the corona radiata than for cortical gray matter, see Fig. 17.
In summary, because of its ultrasoft nature, brain tissue stiffness recordings are highly sensitive to the fluid content of the sample. Undrained samples are stiffer than drained samples, and drainage rates depend critically on the tissue microstructure. These effects are less pronounced in other types of tissues with a lower fluid volume fraction. This explains why the reported stiffness values of brain tissue vary hugely. Without an explicit mention of loading rates, drainage conditions, and sample size and geometry, it is virtually impossible to compare stiffness values recorded under different test conditions. The concept of a single one gray or white matter stiffness value simply does not exist for brain tissue, and it is critical for computational simulations to understand exactly which situation applies to select the appropriate model and parameter values. Figure 6 summarizes various testing techniques to probe the mechanical behavior of brain tissue at different spatial and temporal scales. A prominent method to probe brain tissue at small spatial scales is atomic force microscopy [35]. Unlike nanoindentation on the meso-scale with relatively large indenter tips shown in Fig. 18, atomic force microscopy has resulted in yet higher stiffness in gray than in white matter [94]. In atomic force microscopy, the size of the indenter tip is on the order of the dimensions of individual cells. The indenter seems to be small enough not to trap the fluid beneath the tip which suggests that these tests probe the solid component of brain tissue similar to unconfined experiments at slow loading rates.

Regional Trends Depend on the Length Scale
A prominent method to probe brain tissue at large spatial scales, embedded in the skull, is magnetic resonance Fig. 18 Experimental nanoindentation setup. A 5 mm-thick coronal slice of freshly harvested brain tissue is placed in a 100mm-diameter petri dish and mounted underneath the force transducer of the TriboIndenter TM . A 12 mm-diameter washer marks the indentation region and stabilizes the sample. A circular flat punch of 1.5 mm diameter ensures a homogenized specimen response. Reprinted with permission from [20] elastography [97]. Its major advantage is that it allows us to probe the brain tissue in vivo in its natural environment. Here, we most likely measure additional effects due to the interaction of brain tissue with surrounding structures including the meninges and the skull. Because of its limited spatial resolution, it remains questionable whether in vivo magnetic resonance elastography can accurately quantify regional mechanical properties, especially in thin structures such as the cortex and the corpus callosum [131]. Furthermore, the stiffnesses from magnetic resonance elastography are highly sensitive to the actuation frequency [177]. Not surprisingly, some studies found cortical gray matter to be stiffer than white matter [71], while others reported the opposite trend [114], or no significant differences [10].

Conclusion Towards Regional Trends
In summary, while it seems important to take the significant local variations in tissue properties into account when mechanically modeling brain tissue, we have to pay caution with regard to the specific application we have in mind. There are no right or wrong testing results; depending on the length and time scales during testing, as well as the boundary conditions, we observe different trends. For instance, depending on the application of interest, confined or unconfined testing data could be relevant. The former represent the behavior of brain tissue at intermediate and short time scales during surgery [171] or impact loading [38], while the latter might rather represent the behavior for slow processes during brain development [26,143,157] or tumor growth [90].
In addition to the different regions within the cerebrum, we can distinguish between the cerebrum, cerebellum, corpus callosum, thalamus, and brain stem, amongst others [176]. Indentation experiments, for instance, show that the mouse cerebellum is softer than the mouse cerebrum [108]. This agrees with results on human brain tissue using magnetic resonance elastography [115,123,176]. Finally, it is also important to note that, even within the brain regions we have introduced in Fig. 4, tissue properties may vary significantly. Indentation experiments, for instance, revealed noticeable inter-regional variations within the corona radiata of porcine brain tissue [31] and studies suggest that stiffness variations in white matter tissue are linearly correlated to the local myelin content [173]. Understanding the effects of the tissue microstructure on the macroscopic mechanical response is critical for the interpretation of the constitutive behavior of the human brain for computational simulations.

Open Questions
There is a general agreement that the ultrasoft nature, the high fluid content, and the biochemical composition make brain tissue very different from all other soft tissues in our body. This implies that factors that have traditionally been considered irrelevant in soft tissue mechanics could play an important role when characterizing the material properties of the brain. In this section, we highlight several open questions that could point towards new studies with a view to create a more holistic picture of the behavior of our brain under healthy and diseased conditions.

Is Brain Stiffness Species-Dependent?
Due to limited availability and ethical considerations, only a few studies have actually tested human brain tissue [30,44,50,51,58,59,85,134,151]. Alternatively, researchers consulted porcine [120,121,134,137,139,140,160,165] or bovine brain tissue [13,20,40,163] due to their structural similarities with the human brain. Others tested the properties of rat [35,49,54] or mouse brains [133]. Since the primary goal of developing and calibrating mechanical models for the brain is to assist diagnosis and treatment of human patients, it is important to understand speciesdependent peculiarities. Early compression stress relaxation experiments suggested that monkey brain tissue was stiffer than human brain tissue [59]. More recent indentation experiments show that mouse brain tissue is stiffer than rat brain tissue, which is again stiffer than porcine brain tissue [109]. Interestingly, these observations imply a negative correlation between stiffness and total brain volume: the smaller the brain, the stiffer the response. This hypothesis, however, contradicts shear relaxation experiments in which human brain tissue was stiffer than porcine brain tissue [134]. However, in this study, the specimen thickness of only 1mm might have affected tissue integrity and the obtained results. To date, there is no general agreement on the species-dependence of brain tissue properties and it remains unclear whether the observed differences are an artifact of the testing method or the result of a true size effect that we can observe in microstructural engineering materials.

Is Brain Stiffness Correlated with Cell Density?
One approach to understand and predict regional variations in brain tissue properties discussed in Sect. 3.8 is to disclose the correlation between macroscopic mechanics and the locally varying microstructure. First steps towards this direction have only been taken recently. Due to the functionality of nervous tissue, most attention has been paid to cellular components, while extracellular matrix components were given less consideration. The composition of different cell types such as neurons, astrocytes, oligodendrocytes, or microglia, and even their local morphological changes in response to functional demands is illustrated in Fig. 5. Independent of those differences, different brain regions, and cell decomposition, however, we find a negative correlation between the shear modulus and the total number of cell nuclei, as illustrated in Fig. 19.
These preliminary results agree well with a recent study on live mouse brain tissue [5] and indicate that the cells might actually be the softest component of brain tissue. Accordingly, extracellular matrix components could play an important role. In this context, it is interesting to note that the cells in the central nervous system have been shown to be very soft compared to cells from other tissues [106]-as is the overall tissue response. This conjecture is supported by another recent study showing that the stiffness of brain tissue can not be solely determined by the stiffness of the cells that constitute the tissue [81]. But, it contradicts a previous study on spinal cord tissue, where the stiffness positively correlated with the relative tissue area covered by cell nuclei [94,95]. The latter finding could be attributed to the fact that those measurements were performed using atomic force microscopy indentation on a smaller length scale than the experiments which are the basis for Fig. 19 [22]. On the resolution of individual cells, it was further shown that glial cells, including astrocytes, oligodendrocytes, and microglia, are even softer than neurons [106]. Importantly, however, the stiffness measured for individual cortical cells changed depending on the extracellular matrix material used for coating the dish during experimentation [81]. This emphasizes that when aiming to understand tissue mechanics on the continuum scale, it is insufficient to test the stiffness of each component individually. Rather, the overall tissue response and its correlation with the underlying microstructure needs to be characterized, when cells are embedded in their natural environment.

Is Brain Stiffness Correlated with Myelin Content?
Neurons in the central nervous system are surrounded and cross-linked by myelin, a fatty white substance that wraps around axons to create an electrically insulating layer. Fig. 5 illustrates the microstructural implications of myelin in white matter tissue. While the electrical function of myelin is widely recognized, its mechanical importance remains underestimated. Figure 20 suggests that white matter stiffness is positively correlated with the local myelin content. These results were obtained by combining nanoindentation testing and histological staining in immature, pre-natal brains [172] and mature, post-natal brains [173] and agree well with uniaxial tension experiments on chick embryo spinal cord tissue, which suggest that myelin and cellular coupling of axons via the glial matrix in large part dictate the tensile response of the tissue [150]. The positive correlation between myelination and stiffness is also confirmed by magnetic resonance elastography studies showing that demyelination reduces the stiffness in a murine model of multiple sclerosis [149]. Interestingly, those processes were shown to be reversible after remyelination. We may conclude that myelin is not only important to ensure smooth electrical signal propagation in neurons, but also to protect neurons against physical forces and provide a strong microstructural network that stiffens the white matter tissue as a whole. The strong correlation between the white   [172] matter stiffness and the local myelin content points towards the potential of tissue stiffness as a biomarker for multiple sclerosis and other forms of demyelinating disorders.

Is Brain Stiffness Correlated with DTI Properties?
Microstructural parameters that require histological staining, e.g. the number of cell nuclei and the myelin content, can only be reliably recorded ex vivo. However, it is also desirable to find correlations between mechanical tissue properties and structural data that can be obtained in vivo. One such parameter is the fractional anisotropy (FA) obtained from magnetic resonance imaging and diffusion tensor imaging (DTI). In white matter, fractional anisotropy illustrates the alignment of nerve fibre bundles. In cortical gray matter, it may be interpreted in terms of the orientation and distribution of axonal, dendritic, and glial cell processes [17]. Interestingly, low fractional anisotropy values have recently been attributed to high neuropil complexity [169], which strengthens the hypothesis that brain tissue stiffness closely correlates with interconnections and capillary density [22].  [22]. Importantly, the observed correlation between tissue stiffness and fractional anisotropy could point towards new methods to access regional variations in tissue properties in vivo using magnetic resonance diffusion tensor imaging.

Does Our Brain Stiffen During Development?
Closely linked to the observations in the previous subsections, experimental studies have consistently shown that brain tissue stiffens during development. Our indentation experiments revealed that both bovine gray and white matter tissue stiffened significantly upon maturation: the gray matter stiffness doubled from 0.31 ± 0.20 kPa pre-natally to 0.68 ± 0.20 kPa post-natally; the white matter stiffness tripled from 0.45 ± 0.18 kPa pre-natally to 1.33 ± 0.64kPa post-natally [172]. This is in perfect agreement with a significant increase in the indentation moduli of rat [49,152] and mouse [109] brain tissue beginning at 10-12 days after birth and continuing to 180 days. Interestingly, in exactly this period, myelin basic protein as a measure of the progress of myelination increases, which confirms the close correlation between developmental brain stiffening and myelination discussed in the Sect. 3.9.3. Similarly, dynamic shear experiments on 2-3 day old pig brains yielded significantly lower shear moduli than experiments on adult pig brain samples [160]. Interestingly, also the strain-stiffenig character of the tissue samples increased with maturation [160]. According to magnetic resonance elastography measurements, the adult human brain appears to be three to four times stiffer than the brain of young children [30]. Even adolescent brains still show a softer response than adult brains in certain brain regions including the cerebellum as well as the parietal and temporal lobes [115]. Only one group reported the opposite trend with a decrease in tissue stiffness with age based on indentation and shear experiments on rat brains [62,134].

Does Our Brain Soften with Age?
While it seems well established that brain tissue stiffens during development, the natural next question is whether brain tissue starts to degrade and soften again after it has passed a zenith. Interestingly, neither oscillatory shear tests [51] nor macroscale unconfined experiments [22,85] showed strong age-dependent trends of brain tissue stiffness, as illustrated in Fig. 22. The graphs summarize ex vivo data from ten human brains and indicate that regional trends, as discussed in Sect. 3.8, are markedly stronger than age-or inter-subject-dependent effects-specimens from a specific region yielded moduli in a similar range independent of age or subject [22]. In contrast to these findings, in vivo measurements using magnetic resonance elastography on human subjects indeed yielded a linear decline in whole-brain elasticity within an age range from 18 to 72 years [77,146]. We attribute this observation to changes in the fluid balance of the human brain and a  [22] decrease in total brain volume, which would affect magnetic resonance elastography but not ex vivo measurements on small tissue samples that primarily probe the elastic properties of the solid phase. This hypothesis is supported by the fact that the relative viscous-to-elastic behavior during magnetic resonance elastography did not differ between age groups, suggesting a preservation of the organization of the tissue's microstructure, which is responsible for elastic tissue properties [77].

Does Our Brain Stiffness Change During Disease?
Neurodegenerative diseases involve remodeling of the brain's microstructure. Expectedly, this also leads to changes in the mechanical properties of the tissue. Most insightful in this respect are studies that compare healthy and diseased brain tissue properties in vivo. A limitation of in vivo measurements via magnetic resonance elastography, as already discussed in Sect. 3.9.6, is that they do not necessarily reflect changes in the local stiffness of the solid phase, including cells and extracellular matrix, but rather changes in the overall integrity of brain tissue, including changes in fluid transport and wave propagation. Magnetic resonance elastography shows that brain tissue softens in multiple sclerosis [156,182], Alzheimer's disease [65,124], and demyelination in general [149]. Interestingly, these softening effects scale with disease stage [86]. This points towards an exciting new application of magnetic resonance elastography as a diagnostic tool to diagnose and quantify disease progression. Interestingly, however, while the stiffness seems to be a sensitive marker for taupathology, neuronal loss, and inflammation, this is not the case for amyloid-pathology [65].

Does Our Brain Stiffness Change After Death?
An important unanswered question remains how brain tissue properties measured ex vivo compare to the tissue response in vivo. Several experimental setups have been designed to tackle this issue. Indentation experiments on rats, for instance, showed that the shear modulus obtained in vivo is about 31% higher than that obtained in vitro [118,152]. Similarly, in situ indentation yielded approximately 30 to 50% higher shear moduli than ex vivo indentation [63]. These observations agree well with measurements using ultrasound elastography, where the shear modulus in vivo was about 47% higher than that given by ex vivo measurements [105]. A recent study using magnetic resonance elastography confirms that porcine brain tissue appears stiffer in vivo than ex vivo at frequencies of 100 and 125 Hz [72]. At lower frequencies, however, they found closer agreement between ex vivo and in vivo measurements. Contrary to these finding, other magnetic resonance elastography studies found an increase in shear moduli immediately after death [164], by up to 58% within only three minutes [176]. The origin of this rapid change within such a short period of time is likely of biochemical nature, but has not yet been explored to full extent.

Does Brain Stiffness Change Post Mortem?
Besides death itself, the post mortem storage time could potentially affect experimental results on brain tissue properties. Fig. 23 shows that, when kept intact and hydrated, bovine brain slices maintained their mechanical characteristics from nanoindentation throughout the entire testing period of five days post mortem. Also, in the time window of human brain experiments in Figs. 7, 8, 9, 10 and 11, Fig. 22 Age-dependent shear modulus in different brain regions: cortex (C), basal ganglia (BG), corpus callosum (CC), and corona radiata (CR), and for different loading conditions, simple shear, compression, and tension. Interestingly, independent of age or subject, specimens extracted from the same brain region yield a shear modulus in the same range. All three ex vivo tests display no significant age-dependency we could not observe a notable change in tissue stiffness between samples that were tested first and last. This agrees well with a recent study using ultrasound elastography on Japanese big-ear rabbits, which reported that the change in mechanical properties is negligible at least within 1 hour after death [105]. In contrast, previous studies on porcine brain tissue have revealed a slight increase in tissue stiffness beginning 6h post mortem [61,126]. This may be attributed to the fact that these studies were performed upon continued mechanical loading. During the experiments in Fig. 23, however, we minimized exposure to mechanical testing to clearly separate the effects of mechanical history discussed in Sects. 3.3-3.5 and post-mortem time. We conclude that if the experimental conditions are carefully chosen and the tissue is kept hydrated at all times, the degeneration process of post mortem brains is rather negligible. If brain tissue is stored without any liquid medium, however, the bio-molecular interactions and the mechanical strength of brain tissue deteriorate with prolonged storage duration, for instance due to the degeneration of myelin sheaths and the vacuolization of cristae [183].

Is Brain Stiffness Temperature-Dependent?
Most ex vivo experiments have been performed at room temperature. It is important to understand, how a rise in temperature from room to body temperature in the in vivo situation will affect the mechanical response of brain tissue. A recent study using ultrasound shear wave elastography on rabbit brains indicates that brain tissue stiffness decreases approximately linearly when the temperature increases from room temperature to body temperature, stays relatively constant in the range from 35 to 42 • C , and then rises again [104]. This is an interesting finding as according to those results, the stiffness of brain tissue seems to be constant exactly for the range of temperatures that may occur in vivo. It agrees well with dynamic shear tests on murine brain tissue in which brain tissue was stiffer at 22 • C than at 37 • C [133]. Uniaxial compression experiments on porcine brain tissue, on the contrary, show a slight but insignificant increase in stiffness from room temperature at 22 • C to body temperature at 37 • C [138].

Modeling Aspects
Computational modeling allows us to analyze and predict the behavior of human brain tissue under a variety of loading conditions. However, the value of a numerical prediction critically depends on the adequate choice of constitutive models. In the following, we will systematically propose mathematical formulations to capture the specific characteristics of brain tissue behavior discussed in Sect. 3.
The complexity of the tissue response depends on the loading conditions and so does the appropriate modeling approach. Different constitutive relations may be needed for the same material depending on the particular application. In this section, we will introduce constitutive relations of increasing complexity to capture the elastic, viscoelastic, and poroelastic behavior of brain tissue. We will then make application-specific suggestions towards selecting an appropriate model in the subsequent Sect. 5.
Due to the high compliance of brain tissue and the distinct nonlinearity of the tissue response, even for strains of only 1% as discussed in Sect. 3.1, we limit ourselves to constitutive models using the nonlinear field theory of Fig. 23 Temporal variation of gray and white matter moduli. The consistent moduli within five days post mortem reveal that brain slices are virtually insensitive to the time of preservation. The stiffness increases moderately with indentation depth, from black dots to white dots. Gray matter, left, is consistently softer than white matter, right. Black horizontal lines indicate the mean; gray zones indicate the standard deviation. Adapted from [20] mechanics. To characterize the kinematics of finite deformation, we introduce the deformation map , which maps position X from the undeformed, unloaded configuration, B 0 ∈ ℝ 3 , to its new position x = (X, t) in the deformed, loaded configuration-the current placement of the body at time t, B t [78]. We further introduce the deformation gradient F(X, t) = ∇ X (X, t) to map undeformed line elements to deformed line elements, where X and x denote the position vectors in the unloaded reference and loaded spatial configurations, respectively. The principal stretches a , a = {1, 2, 3} are the square roots of the eigenvalues of the left and right Cauchy-Green tensors defined by b = F ⋅ F t and C = F t ⋅ F , respectively.

Hyperelasticity of Brain Tissue
In a first step, we focus on the time-independent response of brain tissue, neglecting viscous or porous contributions, namely the experimental findings presented in Sects. 3.1 and 3.2 . The main time-independent characteristics are nonlinearity and compression-tension asymmetry. We postulate the existence of a strain-energy function (F) , which is defined per unit reference volume and only depends on the deformation gradient F . We note that several previous studies proposed fiber-reinforced material models for brain tissue, where the strain energy not only depends on the deformation gradient F but also on the fiber direction f 0 [6,37,52,60,66,128,181]. Based on our experimental findings in Sect. 3.7, however, we suggest that the elastic behavior of brain tissue is, to a first approximation, isotropic. Nonetheless, the anisotropy induced by the orientation of nerve fibers may be important for other mechanical processes including damage or diffusion [175].

Hyperelastic Constitutive Modeling
Several phenomenological, isotropic strain-energy functions have been proposed to describe the constitutive behavior of brain tissue [22,41,92,108,137]. Most of these models were originally developed for much stiffer materials such as polymers [130] and calibrated using a single loading mode only [92,108,137]. Here, we evaluate the capability of previously proposed isotropic hyperelastic constitutive models to capture the time-independent response of human brain tissue under multiple loading conditions. We consider three special cases of the generalized Ogden model with the strain-energy function [130], i.e., where the constitutive parameters i correspond to the strainmagnitude-sensitive nonlinear characteristics of the tissue.
The classical shear modulus, known from the linear theory, is given by = 1∕2 ∑ n i=1 i i [78]. Firstly, we adopt the neo-Hookean model with 1 = 2 and 1 = , i.e., Secondly, we use the Mooney-Rivlin model with 1 = 2 , 1 = C 1 = − C 2 , 2 = −2 and 2 = C 2 according to and thirdly we reformulate the one-term Ogden model in terms of the classical shear modulus = 1 1 ∕2 and the parameter = 1 , i.e., In addition, we consider an exponential strain-energy function proposed by Demiray [42] as and a rapidly strain-stiffening material model proposed by Gent [64], Following standard arguments of continuum thermodynamics, we can express the Piola stress tensor P as the derivative of the strain-energy function with respect to the deformation gradient F [78]. Assuming that brain tissue deforms homogeneously and isochorically with the incompressibility constraint det F = 1 , we may provide an analytical prediction of the Piola stresses where n a and N a are the eigenvectors of the left and right Cauchy-Green strain tensors and p serves as a Lagrange multiplier [78]. We can then compare the stresses predicted by the model to experimentally observed responses. Figure 24 illustrates the performance of the hyperelastic constitutive models (2) to (6) to represent the conditioned brain tissue response in different regions, including the cortex, basal ganglia, corona radiata, and corpus callosum during multiple loading conditions, uniaxial compression, uniaxial tension, and simple shear, simultaneously. Table 1 summarizes the resulting region-specific material parameters.

Parameter Identification
As the shear response of brain tissue deviates from linearity, even for small amounts of shear, neither the neo-Hookean nor the Mooney-Rivlin material model are able to satisfactorily represent the experimental data, which is clearly visible in Fig. 24. Only the one-term Ogden model is able to represent all loading modes simultaneously. It not only predicts the nonlinearity but also inherently captures the compression-tension asymmetry with a notably softer response in tension than in compression. These characteristics are controlled by the material parameter : the higher the  Table 1 Constitutive parameters and corresponding coefficients of determination R 2 obtained by calibrating five hyperelastic constitutive models (neo-Hookean, Mooney-Rivlin, Demiray, Gent, and modified one-term Ogden model) with the averaged 'elastic' behavior in differ-ent regions, the corpus callosum (CC), the corona radiata (CR), the basal ganglia (BG), and the cortex (C), for each loading mode separately, and all loading modes simultaneously. Adapted from [22]  absolute value of , the higher the degree of nonlinearity; if > 2 , tensile stresses are stiffer than compressive stresses, if < 2 , we observe the opposite. We note that in simple shear, shear stresses predicted by the constitutive model are independent of the sign of . This implies that we have to be cautious when calibrating the constitutive model exclusively with simple shear data [22]. Figure 25 further demonstrates that the modified oneterm Ogden model is able to capture the inherent characteristic of brain tissue that shear stresses increase significantly with increasing superimposed axial compression but only slightly with increasing axial tension, with coefficients of determination R 2 > 0.91 . This behavior is a logical outcome of the compression-tension asymmetry of brain tissue and can only be captured by one of the strain-energy functions compared in Fig. 24: the modified one-term Ogden model with a negative nonlinearity parameter .
Calibrating the constitutive model with combined compression/tension-shear data yields a similar value for the shear modulus as simultaneously calibrating the model with the data from multiple uniaxial loading modes in Table 1, bottom. However, the absolute value of is much lower. A high absolute value of with ≈ −20 yields unrealistically high shear stresses for high compressive or tensile pre-strain in the combined loading case. In contrast, a low absolute value of with ≈ −7 is not capable of representing the nonlinearity of the shear stress versus amount of shear curve reasonably well. For the sequence of multiple uniaxial loading modes in Fig. 24, the load is limited to 10% strain in compression and tension, and 20% in shear, to not damage the tissue during the course of the experiment. Due to the distinct nonlinearity of the stress-strain curve, even for those relatively small strains, the value obtained from uniaxial experiments would predict unrealistically high stresses for larger strains. This demonstrates that the one-term Ogden model, a phenomenological model in nature, can easily predict an unrealistic behavior when exceeding the deformation range used for parameter identification.
Consistent with these findings, all studies in the literature that considered both compression and tension experiments, reported that only the one or two-term Ogden models could satisfactorily represent the material response [58,121,165]. These studies proposed = 1.0 kPa for cyclic compression-tension experiments on human white matter tissue [58], = 0.8 kPa and = −4.7 for mixed porcine brain tissue [121], and = 0.3-0.7 kPa and = −7.0 when extrapolating tensile porcine white matter data to compression [165]. The lower absolute values for can be attributed to strains of 30% and more.
In contrast, Studies considering each loading mode individually found excellent agreement between experimental data of mixed porcine brain tissue and the Demiray, Gent, and Ogden strain-energy functions [22,137,139,140]. Similarly, a study based on indentation data reported that polynomial, Yeoh, and one-term Ogden models agreed well with experimental data using an inverse finite elements analysis [92]. A more recent study on nanoindentation experiments, in contrast, suggests that the neo-Hookean model best represents indentation data [108].
This emphasizes that, due to the highly complex mechanical response of brain tissue, constitutive models derived from a single loading mode are not necessarily valid for different loading conditions. We conclude that the one-term Ogden model is able to capture the mechanical response of human brain tissue under multiaxial loading modes. However, particular caution is necessary when determining the parameter : The compression-tension asymmetry pre-supposes a negative sign for and high absolute values yield unrealistically high stresses for large strains and multiaxial loading cases.

(In)homogeneous Deformation
The results presented in Figs. 24 and 25, as well as in Table 1, are based on the assumption that brain tissue deforms homogeneously during uniaxial compression and tension as well as in simple shear. However, due to the high compliance of brain tissue and the fact that the upper and lower faces of the specimens are glued to the specimen holders during testing, the deformation actually displays certain inhomogeneities. Figure 26 contrasts the actual deformation of the specimens during simple shear, compression, and tension loadings with finite element simulations using the analytically determined parameters in Table 1. Figure 26 indicates that material parameters calibrated for a homogeneous response tend to underestimate the shear stresses and overestimate the compressive and tensile stresses compared to using inhomogeneous finite element simulations, here exemplary shown for the neo-Hookean and the one-term Ogden model. Especially regarding the oneterm Ogden model, however, the model predictions mostly remain within the standard deviations of the experiments. Figure 25 demonstrates another effect resulting from the inhomogeneous deformation state. Using the analytically calibrated material parameters to numerically simulate combined loading conditions yields lower shear stresses than those predicted analytically, which is in accordance with the results in Fig. 24. However, the simulations capture the qualitative effect that shear stresses not only increase with superimposed compressive strains, but also slightly increase with superimposed tensile strains in accordance with the experimental results in Fig. 25. Taken together, the numerical results in Figs. 25 and 26 emphasize the importance of using an inverse parameter identification scheme to determine appropriate material parameters for computational simulations of brain tissue behavior in the future [167].
In addition to an inverse parameter identification, numerical simulations are valuable to optimize experimental procedures or testing protocols. They allow us-in advance-to evaluate the sensitivity of material parameters towards certain loading conditions, which will help to explicitly design experiments that are suitable for accurate and unambiguous parameter identification.
One effect that can not be captured analytically is the effect of specimen geometry. Figure 27 illustrates how the specimen geometry, height and cross-sectional area, affect the recorded stresses. Stresses significantly increase when specimens become too thin as the deformation inhomogeneities at the fixed faces gain in influence. For compression and tension loadings, a specimen height of approximately 5 mm seems optimal to ensure a consistent response. For simple shear loadings, the effect of specimen height is noticeable independent of the height; however, with increasing height the effect of the cross-sectional area becomes more prominent.

(In)compressibility
The results presented in Fig. 24 and Table 1 are based on the assumption that brain tissue is incompressible, motivated by the high water content of approximately 80%. This assumption may be adequate for impact situations [88,99]; however, especially when considering the time-independent, quasistatic response of brain tissue relevant for extremely slow processes such as brain development, cerebrospinal fluid may escape through the ventricular system, as discussed in Sect. 3.5. This will lead to a slight compressibility, which we can model by adding a volumetric contribution to the modified one-term Ogden model [78,130], where Ogd iso describes the isochoric response, Ogd vol describes the purely volumetric response, ̃a = J −1∕3 a with a = {1, 2, 3} are the volume-preserving parts of the principal stretches, and denotes the bulk modulus.
To demonstrate the influence of compressibility on the brain tissue response during unconfined experiments, we used the material parameters calibrated analytically under the assumption of incompressibility in Table 1 in a finite element setting and varied the Poisson's ratio from 0.3 to 0.49. We ensured that the results for = 0.49 were not affected by locking effects when using linear finite elements through a comparison with the results using mixed finite elements to deal with quasi-incompressibility. Figure 28 shows how tissue compressibility affects the response during unconfined compression, tension, and simple shear. Expectedly, independent of the loading mode, a decrease in the Poisson's ratio also leads to a decrease in tissue stresses. Interestingly, a recent study argues that a different compressibility in compression than in tension leads to the experimentally observed compression-tension asymmetry [167]. From a physical perspective, such an effect could be attributed to poro-elastic effects, which will be discussed in more detail in Sect. 4.3.

Conclusions and Future Perspectives
The one-term Ogden model inherently captures the main characteristics of the time-independent response of brain Fig. 27 Numerical study on the influence of specimen geometry and specimen height h on the recorded tissue response. Stresses show a significant increase for small specimen heights. A rectangular cross-sectional area leads to a direction-dependent response in simple shear, which vanishes when the specimen height approaches zero tissue, the nonlinearity and compression-tension asymmetry. It is capable of representing multiple loading modes simultaneously. However, material parameters calibrated analytically assuming incompressibility and a homogeneous deformation will tend to overestimate compressive and tensile stresses and underestimate shear stresses when used in finite element simulations. To address these limitations, we recommend calibrating the model parameters with sophisticated inverse identification schemes, which capture inhomogeneous deformation states and several loading modes simultaneously. A remaining open challenge is to identify a single model that captures a wide range of strains. Depending on the loading conditions, it may be appropriate to model human brain tissue as an incompressible or compressible solid [88,99]. If movement of free flowing cerebrospinal fluid into the ventricles and the subarachnoid space is possible, for instance during slow processes such as brain development, brain tissue effectively changes its local volume [148]. Therefore a slight compressibility with a Poisson's ratio of 0.45 to 0.49 seems appropriate [26]. In addition, parameters should be calibrated using the conditioned, drained tissue response. In impact situations, in contrast, the fluid offers resistance and contributes to tissue stiffness. In these case, the tissue may behave incompressibly and parameters should be calibrated using the unconditioned tissue response.

Finite Viscoelasticity
While hyperelasticity can describe the nonlinearity and compression-tension asymmetry as outlined in Sects. 3.1 and 3.2 , it cannot capture the characteristic time-dependent effects, the highly hysteretic behavior of Sect. 3.3 [22,58,121], the significant strain-rate-dependence of Sect. 3.4, and the characteristic conditioning as pointed out in Sect. 3.5. The time-dependent nature of brain tissue is associated with various physical mechanisms and time scales: The motion of fluid within the solid network of cells and extracellular matrix introduces a poroelastic behavior [48,58], whereas intracellular interactions between the cytoplasm, nucleus, and the cell membrane trigger a viscoelastic response [84]. Importantly, when tested below a critical threshold, the conditioning behavior of brain tissue is fully recoverable, as we conclude from Fig. 15: No damage or permanent softening occurs, which justifies to model the time-dependent response of brain tissue using visco-or poro-viscoelasticity instead of irreversible damage or the Mullin's effect [58].
In this section, we will extend the hyperelastic model presented in Sect. 4.1 to account for time-dependent, viscoelastic effects. The simplest and most popular approach to characterize the time-dependent behavior of brain tissue is to use a Prony series [20,57,100,120,121,134,139]. A Prony series approach is equivalent to a generalized Maxwell model for linear viscoelasticity in relaxation type loading [89]. However, the Prony series approach has two major limitations: it is restricted to linear elasticity and is computationally expensive [41]. To account for the large deformability of brain tissue discussed in Sect. 2.1, instead, we consider a class of viscoelastic models within the general setting of finite deformation continuum mechanics. We multiplicatively decompose the deformation gradient into elastic and inelastic parts [153], additively decompose the free-energy function into an equilibrium and non-equilibrium parts [154], and introduce internal variables to account for the rate-dependent behavior [142].

Kinematics of Finite Viscoelasticity
To model the viscoelastic nature of brain tissue, we decompose the deformation gradient F into elastic and viscous parts, where i denotes the parallel arrangement of m viscoelastic elements [153], as exemplary shown for two viscoelastic  Fig. 29. We can then introduce the spatial velocity gradient as and decompose it additively into elastic parts, l e i =Ḟ e i ⋅ (F e i ) −1 , and viscous parts, From the deformation gradient, we determine the left Cauchy-Green tensor b and its spectral representation in terms of the principal stretches a and the Eulerian eigenvectors n a , It is convenient to introduce the elastic left Cauchy-Green tensor for each mode, i.e. with eigenvalues ( e i a ) 2 and eigenvectors n e i a , which are, in general, not identical to the eigenvectors of the total left Cauchy-Green tensor, n e i a ≠ n a . The material time derivative of the elastic left Cauchy-Green tensor b e i is which introduces its Lie-derivative, i.e.
2 a n a ⊗ n a .
along the velocity field of the material motion.

Viscoelastic Constitutive Modeling
Motivated by the experimental findings discussed in Sect. 3.7, we assume an isotropic, incompressible material response for both the elastic and the viscoelastic behavior. We introduce the viscoelastic free-energy function as the sum of three terms, an equilibrium part eq in terms of the total principal stretches a , a non-equilibrium part neq = ∑ m i=1 i in terms of the i = 1, … , m elastic principal stretches e ia , and a term p [ J − 1 ] that enforces the incompressibility constraint J − 1 = 0 via the Lagrange multiplier p [78], i.e.
Similarly, we introduce the stress power P as the sum of an equilibrium part P eq = eq ∶ l in terms of the equilibrium Kirchoff stress eq and a non-equilibrium part P neq = neq ∶ l in terms of the non-equilibrium Kirchoff We can then evaluate the dissipation inequality , in terms of the Lie derivative of the elastic left Cauchy-Green tensor (14), and obtain the following explicit representation of the dissipation inequality The first term of (18) vanishes when we define the equilibrium Kirchhoff stress tensor as (17) D = [ P eq + P neq ] − [̇e q +̇n eq ] ≥ 0.
eq a a n a ⊗ n a , Fig. 29 Multiplicative decomposition model, where F is associated with the main elastic network characterized by material parameters ∞ and ∞ , F v 1 is the viscous damper associated with fluid flow inside the cell characterized by the viscosity 1 with the corresponding hyperelastic spring F e 1 with parameters 1 and 1 , and F v 2 is the viscous damper associated with fluid flow within the solid network of cells and extracellular matrix characterized by the viscosity 2 with the corresponding hyperelastic spring F e 2 with parameters 2 and 2 . Adapted from [23] and the second term vanishes for non-equilibrium Kirchhoff stresses, i.e. This leaves the reduced dissipation inequalities for each individual mode i, which is It remains to specify the equilibrium and non-equilibrium parts of the free energy, eq and neq , and the evolution of the internal variables b e i . For the equilibrium energy, we adopt a modified one-term Ogden model parameterized in terms of the total stretches a [130], which best represented the hyperelastic response of human brain tissue in Sect. 4.1, i.e.
with the equilibrium shear modulus ∞ and the equilibrium compression-tension asymmetry parameter ∞ [22]. The derivative in (19)

then becomes
For the non-equilibrium energy, we adopt the same Ogden strain-energy function which introduces two additional parameters i and i for each mode i, and is now parameterized in terms of the deviatoric i a n e i a ⊗ n e i a .  {1, 2, 3} and a ≠ b , a ≠ c , and b ≠ c . Finally, we need to specify the temporal evolution of the viscoelastic kinematics. To a priori satisfy the reduced dissipation inequality (21) [142,154,170], and motivated by the linear relation between hysteresis and maximum stress during our experiments, as illustrated in Fig. 30, we choose the evolution equations for the internal variables b e i as linear functions of the Kirchhoff stresses i [18,79], i.e. This introduces one additional parameter for each mode i, the viscosity i > 0 , or, when scaled with the corresponding shear modulus i , the associated relaxation time, i = i ∕ i [23,142]. Since the internal variables b e i are a linear function of the deviatoric Kirchhoff stress tensors i , the elastic deformation always remains volume preserving, J e i = 1 , and thus J v i = J∕J e i = 1 . To advance the internal variables in time [70,80,142], we may adopt an implicit or explicit time integration scheme [23]. Finally, we calculate the Piola stress P as the right-sided pull back of the Kirchhoff stress , i.e.
We note that we may decompose the viscous parts of the spatial velocity gradient,

into their symmetric and skew-symmetric parts, l
, and the viscous spin rates, If we now reformulate

Fig. 30
Hysteresis versus maximum stress during the first and third loading cycles for all specimens tested in shear, compression, and tension. The data indicate that the dissipated energy increases linearly with the maximum recorded stress. Adapted from [23] the reduced dissipation inequalities (21) using the symmetry of i , we obtain A comparison with equation (21) yields an alternative expression for the evolution equations for the internal variables (26), expressed in terms of the viscous stretch rates which has been proposed by [80]. If all the rheological elements in the model were linear, it would reduce to the generalized Maxwell model, which is equivalent to the Prony series in a classical relaxation experiment [41,155]. Figure 31 demonstrates that the proposed viscoelastic constitutive model is capable of capturing the experimentally observed pronounced pre-conditioning during the first loading cycle and the moderate conditioning effects during all subsequent cycles [23]. The model also predicts the experimentally observed successive softening when increasing the strain in a step-wise fashion [23], as discussed in Sect. 3.5. In simple shear, however, the model predictions deviate from the experimental data, especially for the initial loading segment. Furthermore, the model predicts larger residual stresses upon completing the first loading cycle of each stretch level than observed in the experiment. This difference is likely caused by the pore fluid that squeezes out of the sample when it is first loaded [58]. Our monophasic viscoelastic model can only implicitly capture these porous effects [11]. To accurately model the fluid flow within the tissue, we will introduce a biphasic poro-viscoelastic model [33,48,68] in the next section to further contrast and compare viscous and porous effects in human brain tissue. Figure 32 demonstrates the effect of the strain rate on the observed stress response. The viscoelastic model nicely captures the experimentally observed increase in the maximum stresses with increasing strain rate, as discussed in Sect

Parameter Identification
We have shown that the finite viscoelastic model can qualitatively capture the most important time-dependent characteristics of brain tissue including hysteresis, Fig. 31 Representative conditioning behavior for cyclic loading with step-wise increasing maximum strain. Shear/nominal stress versus amount of shear/stretch behavior for specimens from the corona radiata; experiment (top row) and model prediction (bottom row). At a stretch level of 0.9, the difference between the first and second cycles indicates pronounced pre-conditioning. When increasing the stretch level to 0.8, the first curve initially follows the pre-conditioned behavior up until 0.9; beyond 0.9, the curve again displays pronounced pre-conditioning between the first and second cycles. The model is capable of predicting this history-dependence at different load levels. For simple shear loading, the qualitative model prediction slightly deviates from the experimental observations strain-rate-dependence, and conditioning. In order to use the proposed model in finite element simulations, we need to provide material parameters that also quantitatively capture the tissue response. Instead of following the common approach of previous studies [57,120,137] and calibrating the elastic parameters using only the initial loading path during quasi-static experiments and the viscous parameters using only the stress relaxation behavior, we integrate data from simple shear, unconfined compression, tension, shear relaxation, and compression relaxation experiments [22] to holistically characterize the overall behavior of brain tissue [23]. Notably, all tests were performed on one and the same samples to avoid that inter-specimen variations affect the results. By simultaneously considering the data of all five tests, we identified material parameters that characterize the behavior of the exact same tissue samples under all five loading modes. We further use the entire loading history for parameter identification-including unloading during cyclic experiments and the loading ramp for stress relaxation experiments [23]. Figures 33 and 34 demonstrate that the finite viscoelastic constitutive model with two viscoelastic elements is capable of representing the average experimental data during the first and third loading cycles associated with the unconditioned and conditioned response, respectively, for the four different brain regions, the cortex, the basal ganglia, the corona radiata, and the corpus callosum [23]. Tables 2 and  3 summarize the corresponding region-specific material parameters. The model slightly underestimates the shear response, especially during the initial loading segment, and slightly overestimates the stresses under tensile loading but nicely captures the relaxation behavior over the entire loading history. Figure 35 contrasts the shear moduli of the unconditioned and conditioned tissue responses in four brain regions [23]. Table 4 summarizes the corresponding characteristic time constants i = i ∕ i near thermodynamic equilibrium. Shear moduli are generally lower for the conditioned than for the unconditioned response, which may be attributed to the contribution of free flowing fluid, which only offers resistance during initial loading but not during subsequent cycles.
The equilibrium shear modulus ∞ displays similar regional trends for the unconditioned and conditioned responses. It is highest in the cortex, lowest in the corpus callosum, and comparable in basal ganglia and corpus callosum. Interestingly, we observe the same regional trends for the shear modulus 2 of the conditioned response in the Maxwell element with particularly high time constants 2 . In contrast, the shear moduli 1 corresponding to lower time constants 1 are generally higher than ∞ and 2 , and show slightly different regional trends: White matter regions, corona radiata and corpus callosum, have stiffened relative to gray matter regions, cortex and basal ganglia.
It is interesting to note that the conditioned response displays higher time constants 2 than the unconditioned response, although we have calibrated both parameter sets with exactly the same stress relaxation experiments. This indicates that the cyclic experiments have a pronounced effect on the viscoelastic parameter identification. Comparing the unconditioned and conditioned behavior suggests that we can attribute the lower time scale 1 , to the viscous component of the solid phase, and the higher time scale 2 , to porous effects of the fluid phase: At short time scales, fluid moves inside the cell; at intermediate time scales, fluid moves through the porous solid skeleton of cells and extracellular matrix, as sketched in Fig. 29 [19]. When using the conditioned response of the third cycle for our calibration, we intentionally neglect the porous effect and, accordingly, the higher time constant adopts significantly larger values than for the unconditioned response. This agrees well with a previous study that has reported a pre-conditioned viscosity of 60 kPas in unconfined compression tests of porcine brain [135].

Conclusions and Future Perspectives
We have proposed a finite viscoelastic model that combines the hyperelastic Ogden model with two viscoelastic elements and can, in addition to the experimentally observed compression-tension asymmetry and nonlinearity, capture timedependent effects including hysteresis according to Sect. 3.3, and the successive softening for stepwise loading according Fig. 33 Simultaneous identification of viscoelastic material parameters for the finite viscoelastic constitutive model with two generalized Maxwell elements and a single nonlinearity parameter = i = ∞ for all five loading conditions. Average experimental data (solid lines) during the first loading cycle associated with the unconditioned response with standard deviations in four regions, cortex, basal ganglia, corona radiata, and corpus callosum with corresponding constitutive model (dashed lines) calibrated using data from all loading modes simultaneously, see Table 2. Adapted from [25] to Sect. 3.5. Importantly, we have shown that we can provide a single set of material parameters that simultaneously fits cyclic loadings and the early and late relaxation behavior of brain tissue. We have provided separate parameter sets for the unconditioned response using the first loading cycle and for the conditioned response using the third loading cycle. Expectedly, the unconditioned tissue response is markedly stiffer than the conditioned response. We would like to Fig. 34 Simultaneous identification of viscoelastic material parameters for the finite viscoelastic constitutive model with two generalized Maxwell elements and a single nonlinearity parameter = i = ∞ for all five loading conditions. Average experimental data (solid lines) during the third loading cycle associated with the conditioned response with standard deviations in four regions, cortex, basal ganglia, corona radiata, and corpus callosum with corresponding constitutive model (dashed lines) calibrated using data from all loading modes simultaneously, see Table 3. Adapted from [25] point out that there is no right or wrong set of parameters: Depending on the application of interest-for example the interpretation of an ex vivo test or the prediction of an in vivo response-either the unconditioned or conditioned parameter set seems to be more physiological.
Through our viscoelastic parameter identification, we have obtained insights that help us to better understand the individual time-dependent contributions of the solid and fluid phases [58]. The first viscoelastic element corresponds to time constants on the order of 1 to 3 seconds and can be associated with the viscoelasticity of individual cells, while the second viscoelastic element corresponds to time constants between 20 and 100 seconds and can represent the movement of fluid within the solid skeleton of cells and extracellular matrix. Such insights become especially important when we attempt to interpret and model brain tissue as a biphasic material [33,48,57], which will be the subject of the next section.
While our simultaneous analysis of multiple loading conditions has significantly improved the parameterization compared to previous viscoelastic models based on a single loading mode [13,80,135], the asymmetry predicted by our model is still less pronounced than observed in our experiments. Furthermore, the parameters are based on the assumption of a homogeneous deformation state. To improve the fit of the model, an inverse parameter identification scheme could be used in the future to potentially address these concerns. Further open challenges regarding the viscoelastic modeling of brain tissue remain the discrepancies between model and experiment under shear loading and the overestimation of the tensile stresses. Table 2 Viscoelastic parameters and coefficients of determination for the viscoelastic constitutive model with a single nonlinearity parameter = ∞ = i calibrated with the averaged unconditioned response in four brain regions, the cortex (C), the basal ganglia (BG), the corona radiata (CR), and the corpus callosum (CC), for different loading modes simultaneously, simple shear, compression, tension, shear relaxation, and compression relaxation, see Fig. 33. Adapted from [25] unconditioned  Table 3 Viscoelastic parameters and coefficients of determination for the viscoelastic constitutive model with a single nonlinearity parameter = ∞ = i calibrated with the averaged conditioned response in four brain regions, the cortex (C), the basal ganglia (BG), the corona radiata (CR), and the corpus callosum (CC), for different loading modes simultaneously, simple shear, compression, tension, shear relaxation, and compression relaxation, see Fig. 34. Adapted from [25] conditioned

Finite Poro-viscoelasticity
To explicitly model the biphasic nature of brain tissue, we can combine the finite viscoelastic model for brain tissue introduced in Sect. 4.2 with the theory of porous media [47]. Importantly, a poro-elastic model alone would not be sufficient to model the hysteretic behavior of brain tissue [39].

Kinematics of Poro-viscoelasticity
We assimilate brain tissue to a biphasic material consisting of an immiscible aggregate of a solid skeleton (S) fully saturated with pore fluid (F). Both components are assumed to be incompressible. The solid phase represents the network of cells and extracellular matrix-including the fluid trapped inside the cells and within the extracellular matrix-while the pore fluid exclusively represents the free-flowing fluid in brain tissue. Figure 36 illustrates that, at a given time t, the constituents simultaneously occupy the same spatial position x in the current configuration B t . The primary unknowns of finite poro-viscoelasticity are the solid deformation x = (X S , t) and the fluid pore pressure p. The strong form of the governing equations in the reference configuration reads where the Piola stress P(F S , p) and the seepage velocity in the reference configuration W(F S , p) are functions of the solid deformation gradient F S = x∕ X S , and the fluid pore pressure p, B 0 is the body force vector in the reference configuration, and J S = detF S is the solid Jacobian. The seepage velocity describes the velocity of the fluid phase with respect to the deforming solid. Following the considerations in the previous section, we involve finite viscoelasticity by decomposing the solid deformation gradient into elastic and viscous parts, F S = F e S ⋅ F v S , and introduce the elastic left Cauchy-Green tensor b e = F e S ⋅ (F e S ) t = ∑ 3 a=1 ( e a ) 2 n e a ⊗ n e a , with eigenvalues ( e a ) 2 and eigenvectors n e a , a = 1, 2, 3.

Constitutive Modeling-Solid Component
We further additively split the strain-energy function of the solid skeleton into an equilibrium part eq and a non-equilibrium part neq , S = eq + neq . The total free-energy function then reads = eq ( i ) + neq ( e i ) + U(J S ) , where U(J S ) is an extension function that describes the compressibility effects of the poroelastic material, including the concept of a compaction point. We select a function that introduces the Lamé parameter of the solid component, where n S0 is the initial porosity of the biphasic material, i.e., the volume fraction of the solid component in the reference

Constitutive Modeling-Fluid Component
In addition to the constitutive equations for the solid, we need to specify a material model for the fluid defined by the seepage velocity W = w ⋅ cof F . Standard approaches have used Darcy's law to describe the motion of the pore fluid with respect to the deforming solid component in terms of the changing pressure gradients in the biphasic material and its deformation-dependent permeability [47,113]. In this case, the seepage velocity is defined as where F and F are the effective fluid viscosity and fluid density, K S is the intrinsic permeability tensor, and g is the vector of gravity acceleration. For example, we can assume an isotropic, deformation-dependent intrinsic permeability tensor K S = ([J S − n S0 ]∕[1 − n S0 ]) K S0 , where the initial intrinsic permeability tensor K S0 = K S0 I would reduce to the weighted identity tensor when the material permeability is Fig. 36 Kinematics of the material configuration B 0 , and the spatial configurations B t and B of a biphasic material body within the context of finite poroelasticity [39] isotropic [113]. Unfortunately, our current experimental data are not sufficient to adequately calibrate the poro-viscoelastic constitutive model. To truly calibrate a poro-viscoelastic model, we would have to perform experiments with drained and undrained conditions [58] to quantify the-potentially region-and loading-mode-dependent-permeability of brain tissue and to clearly differentiate between viscous and porous effects.

Conclusions and Future Perspectives
While the finite viscoelastic constitutive model introduced in Sect. 4.2 can implicitly capture poro-elastic effects through its two independent time scales, a poro-elastic model alone is insufficient to describe the hysteretic behavior of brain tissue [99]. The advantage of using a more complex biphasic poro-viscoelastic model is that it allows us to explicitly control the physics of the fluid as it moves within the solid skeleton of cells and extracellular matrix. Notably, we could even further advance the poro-viscoelastic constitutive model by adding components such as the blood, which flows inside the small brain capillaries [48]. In this case, brain tissue is permeated by two viscous pore liquids, the interstitial fluid and the blood. An important drawback concerning the poro-viscoelasticity of brain tissue is that, to date, the poro-elasticity of brain tissue has not yet been well characterized experimentally. Consolidation tests performed on human white matter brain tissue have shown that brain tissue behaves similiarly to filled elastomers [58]. However, the data available in the literature are not sufficient to adequately calibrate the poro-viscoelastic model presented in Sect. 4.3. Well designed experiments will be essential to quantify region-specific permeability of brain tissue and to differentiate between viscous and porous effects. This will be critical when choosing the appropriate constitutive model for a different biomedical application.

Application-Specific Considerations
Computational modeling can form a cornerstone in improving diagnostics and treatment of brain tissue injury and disease. However, reliable, clinically relevant predictions of brain tissue behavior critically depend on the appropriate choice of constitutive models and material parameters for the numerical simulations. Importantly, the required model complexity depends on the application of interest and the corresponding length-and time-scales. Brain deforms differently during trauma at high rates, surgery at medium rates, and brain development at low rates. To save computational cost, we should evaluate efficiency versus accuracy, and choose a model that is as complex as necessary but as simple as possible. In the following, we present selected examples in which numerical simulations can valuably assist medical diagnostics and treatment of neurological disorders. We discuss the relevant length-and time-scales and propose an appropriate modeling approach. We intend to demonstrate how the evaluation of relevant scales, as summarized in Fig. 6, can help to choose the appropriate mechanical model for a certain application.

Brain Development
Cortical folding is a classical hallmark of most mammalian brains and closely correlates with brain function and dysfunction. During development, the initially smooth surface evolves into an elaborately convoluted pattern. Growing evidence suggests that mechanical instabilities emerging from differential growth between a faster growing outer gray matter, the cortex, and a slower growing inner white matter, the subcortex, play a major role in brain morphogenesis. We can model brain growth by complementing the kinematics of finite deformation introduced in Sect. 4 by the kinematics of finite growth [144]. We multiplicatively decompose the deformation gradient F into an elastic part F e and a growth part F g [144], F = ∇ X = F e ⋅ F g . A key challenge is to prescribe the growth tensor F g and its evolution in time to realistically mimic the underlying biological processes [3] in the developing brain. In the simplest case, growth is purely morphogenetic, evolving in time irrespective of the mechanical environment [117]. In more complex cases, growth is mechanically driven, for instance by hypertension in arteries [98] or by hyperstretch in skin [159].
Inspired by cellular mechanisms during brain development [27], we assume that both cortex and subcortex grow isotropically with F g = I , where is a scalar-valued growth multiplier. The grown volume J g = det F g = 3 is identical to the growth multiplier cubed 3 . We further assume that cortical growth is purely morphogenetic [2]independent of mechanical stress or strain [4]-growing linearly in time at the growth rate G c = const, and that subcortical growth is stretch-dependent at a growth rate G s , as illustrated in Fig. 37 [21,26]. This yields the evolution of the cortical growth multiplier, ̇c = G c and the subcortical growth multiplier,

The term in the Macaulay brackets
⟨ J e − J 0 ⟩ activates growth only if the elastic volume stretch J e exceeds the threshold value J 0 , i.e., when axons are stretched beyond their physiological limit [43]. Taken together, the mechanical model of brain growth enables us to bridge the scales from accumulation and growth of individual cells into the growth kinematics of cortex and subcortex, and ultimately cortical folding [67]. Figure 38 demonstrates that the model explains why the surface-to-volume ratio in different mammalian brains increases disproportionally with brain size and why brains of larger mammals tend to be more folded than brains of smaller animals [179]. From a mechanics point of view, this phenomenon is a merely geometric effect and independent of phylogenetic effects-contrary to the common notion in developmental biology [179]. The results in Fig. 38 demonstrate the great opportunities of using mechanical models to understand brain development. However, the emerging folding pattern is the result of an instability problem, which, by its very nature, is highly sensitive to the constitutive model and material parameters of finite element simulations.
Concerning the choice of the appropriate modeling approach, brain folding occurs over weeks and months of gestation and we can safely assume that the cerebrospinal fluid has time to escape through the ventricular system and that its contributions to the tissue stiffness do not change markedly during the simulation period. Consequently, we propose to use an isotropic hyperelastic constitutive model calibrated with the conditioned experimental response. However, we need to pay caution with regard to the material parameters since extremely large deformations will occur in the bottom of sulci-the valleys of cortical folds. If we use the material parameters calibrated using experimental data up to 10% compressive and tensile strains, and 20% shear strains from Table 1, the model will predict unrealistically high stresses in the valleys of emerging folds and unrealistic folding patterns as evident in the bottom row of Fig. 39. This is a natural result of the corresponding nonlinearity parameter = −20 . Even for a slightly reduced material nonlinearity with = −7 in the center row, according to the parameters identified using combined compression/tensionshear experiments, the folding pattern still deviates from the one observed in the real brain. Surprisingly, we obtain the most realistic patterns for the neo-Hookean material model with = 2 in the top row of Fig. 39. These results emphasize that we need to be careful when using material parameters for deformation states that exceed the experimentally investigated ranges of strains. For the application of cortical Continuum model for cortical and subcortical growth. The cortex, the gray matter, grows morphogenetically at a constant rate G c . Cortical growth induces subcortical deformation, which triggers subcortical growth. The subcortex, the white matter, grows at a stretch-dependent rate as G s ⟨ J e − J 0 ⟩ , where G s mimics the axon elongation rate and ⟨ J e − J 0 ⟩ activates growth only, if the elastic volume stretch J e exceeds its baseline value J 0 . Adapted from [28] Fig. 38 Variety of mammalian brains. Mammalian brains hugely vary in size, shape, and convolutional complexity but only little in cortical thickness (top row), adapted from [179]. Numerical simulation of a growing shell on a spheroidal substrate for increasing ratio of elliptical radius to cortical thickness R / t (bottom row). Folding complexity increases with absolute size Fig. 39 Sensitivity of the evolving surface morphology during brain development towards the constitutive model used for brain tissue. With increasing absolute value of the nonlinearity parameter from top to bottom, as introduced in Sect. 4.1, the wavelength increases and the primary pattern transitions from sinusoidal wrinkles to cusped sulci folding, where high compression occurs locally in the developing sulci, the neo-Hookean model indeed seems to be the best choice. Alternatively, in a more physiological model, we could account for the continuous adaption of the living brain tissue to mechanical loading-for example, motivated by the accumulated cell death in the bottom of the sulciwhich gradually reduces the tissue stresses in these regions.

Neurodegenerative Diseases
Because of demographic changes worldwide, neurodegeneration will, without a doubt, become a major challenge in medicine and public health. More than 45 million people are living with dementia today and this number is expected to triple by 2050 [1]. A major challenge of neurodegenerative diseases is that the initial symptoms of cognitive decline become noticeable only one or two decades after the first pathological abnormalities have started to affect the brain [166]. For several decades, researchers have tried to establish early predictive indicators of neurodegeneration and identify common pathological themes of neurodegenerative disorders. Computational simulations are now increasingly recognized as a powerful tool to correlate biochemical and biomechanical events and provide early biomarkers for neurodegenerative disorders [175]. Figure 40 shows how we can use a physics-based reactiondiffusion model to identify common underlying schemes in the pathogenesis of Alzheimer's disease, Parkinson's disease, and amyotrophic lateral sclerosis [174]. The underlying simulation combines the classical Fisher-Kolmogorov equation [55], ̇c = Div (D ⋅ ∇c) + c [ 1 − c ] for the toxic protein concentration c with an anisotropic diffusion model, D = d ext I + d axn n ⊗ n , with a pronounced diffusion along the axonal direction n . Here is the growth rate of toxic protein, and d ext and d axn are the diffusion parameters associated with extracellular diffusion and axonal transport, respectively. The model correctly predicts amyloid-deposits and tau inclusions in Alzheimer's disease, -synuclein inclusions in Parkinson's disease, and TDP-43 inclusions in amyotrophic lateral sclerosis and displays excellent agreement with the histological patterns observed in diseased human brains [87]. When integrated across the brain, the concentration profiles result in biomarker curves that display a striking similarity with the sigmoid shape and qualitative timeline of clinical biomarker models [82]. Regions of high toxic protein concentration in Fig. 40 are likely correlated with neuronal death and tissue atrophy, which we can model as negative tissue growth, F = F e ⋅ F a . To model the biochemical effects of neurodegeneration, we could introduce the atrophy tensor F a as a function of the toxic protein concentration c [175]. To model the biomechanical effects of neurodegeneration, we recommend using a finite hyperelastic model and assume that viscoelastic and poro-elastic effects are negligible on the relevant time scales on the order of several decades [75]. While the computational modeling of neurodegeneration is still in its infancy, there seems to be a general agreement that a more quantitative understanding of the spatio-temporal spreading of neurodegenerative diseases is necessary to establish a prognostic timeframe of disease progression. This could have important clinical implications, ranging from more accurate estimates of the socioeconomic burden of neurodegeneration to a more informed design of clinical trials and pharmacological intervention.

Hydrocephalus
Hydrocephalus is a medical condition in which excess cerebrospinal fluid builds up in the ventricles of the brain. The excess fluid can increase the size of the ventricles and, with it, the intracranial pressure. Normal pressure hydrocephalus develops gradually in time and manifests itself in enlarged ventricles, while the intracranial pressure remains at normal physiological levels [73]. The common treatment strategy is cerebrospinal fluid flow diversion and surgical placement of a shunt to drain excess fluid into the abdomen. The improvement after such a neurosurgical shunt procedure varies greatly with a morbidity of about 40-50% Fig. 40 Spatio-temporal progression of toxic protein across the brain for different initial seeding regions associated with a, b Alzheimer's disease, c Parkinson's disease, and d amyotrophic lateral sclerosis. Adapted from [175] [102]. Therefore, there is an urgent need for numerical simulations that assist the design of better diagnostic and treatment protocols.
Efforts towards this direction have employed mostly twodimensional finite element models using poro-elasticity with a linear elastic solid skeleton [88,103,125], or the multiplenetwork poroelastic theory [161]. Regarding the conclusions in Sect. 4.3, we infer that it may be important to also take viscoelastic effects into account, which has been done in [34]. Most recent works even propose a triphasic mechano-electrochemical theory [45]. This model additionally accounts for a change in the ionic concentrations of the ventricular cerebrospinal fluid in the absence of an elevated intracranial pressure as a possible cause for the onset of normal pressure hydrocephalus. Interestingly, a three-dimensional finite element study found no major advantage by using a biphasic instead of a single phase model for brain tissue [46].
Based on the different time scales-rather rapid in hydrocephalus in young children versus gradual in normal pressure hydrocephalus in the elderly-either a poro-viscoelastic or a purely poro-elastic model could be more appropriate versus more efficient. It remains to carefully evaluate the influence of viscous and porous effects to decide which model to use. Intracranial pressure monitoring is critical in many neurological conditions. The current gold standard is to monitor the intracranial pressure via catheter transducer systems that are inserted into the cranium. This process is highly invasive and there is hope that computational simulations can help estimate the intracranial pressure less invasively.

Drug Delivery
A promising treatment approach for neurological disorders is the infusion of therapeutic agents, generally known as convection-enhanced drug delivery [16]. Infusion catheters are placed through small holes in the skull directly into the extra-vascular space of the brain parenchyma. The pressure gradient generated by external medical pumps initiates an interstitial fluid flow and, therewith, the distribution of the therapeutic agents. Various poro-elastic mechanical models have been developed to simulate this problem analytically [122] or numerically [32,48,155]. Drug delivery is a typical example where it seems indispensible to model the fluid phase of brain tissue explicitly by using a biphasic constitutive model based on the findings in Sect. 4.3 and due to the fact that infusion of therapeutic agents occurs at intermediate time scales, it might also be important to include viscous effects by using a poro-viscoelastic model.

Neurosurgery
During neurosurgical procedures, brain tissue can temporarily experience strains and stresses that are high enough to induce tissue damage. In this respect, simulation tools capable of predicting the level of stress and strain can be highly valuable to optimize surgical procedures [99]. A typical example is decompressive craniectomy, a traditional but controversial surgical procedure that removes part of the skull to allow an injured and swollen brain to expand outward and release an elevated intracranial pressure [93]. Despite a marked increase in popularity, a decompressive craniectomy remains a compromise between maximizing control of the intracranial pressure and minimizing side effects to the surrounding tissue [69]. To date, the precise criteria when to perform a decompressive craniectomy, the optimal timing of treatment, the optimal location and size of the skull opening, and the long-term functional outcome remain unclear [171]. Studies suggest that mechanical strains beyond the physiological level are associated with its undesired, high failure rates. However, the precise strain fields induced by the craniectomy are unknown. Figure 41 shows finite element simulations of a decompressive craniectomy performed on a personalized head model from magnetic resonance images [178]. The simulations compare the two most common types of craniectomy, a unilateral flap and a frontal flap. For each case, either both hemispheres or the left or right hemisphere of the brain are swollen and part of the skull is removed to allow the brain to bulge outward and release the intracranical swelling pressure. Clearly, the underlying swelling is a biochemical process associated with the absorption of fluid, and we can best model this phenomenon using a poro-viscoelastic approach. While swelling takes place on the order of hours, the clinical procedure of skull opening takes place on the order of minutes on which both porous and viscous effects can become relevant. The simulations in Fig. 41 reveal three potential failure mechanisms associated with the procedure: axonal stretch in the center of the bulge, axonal compression at the edge of the opening, and axonal shear around the opening [171]. Strikingly, for a swelling of only 10%, axonal strain, compression, and shear reach local maxima of up to 30%, and exceed the reported functional and morphological damage thresholds on the order of 20% [8]. Computational simulations of surgical procedures can help quantify brain deformation, tissue strain, axonal stretch, and shear with the goal to identify high-risk regions for brain damage on a personalized basis [56]. While computational modeling is beyond clinical practice in neurosurgery today, simulations of neurosurgical procedures have the potential to rationalize surgical process parameters including timing, location, and size, and provide standardized guidelines for clinical decision making and neurosurgical planning.

Traumatic Brain Injury
A major public health concern are the consequences of neurotrauma with over 2 million people affected by traumatic brain injury each year [111]. During traumatic brain injury, external mechanical load leads to damage of the highly delicate brain tissue and ultimately loss of brain function [116]. Clinically, traumatic brain injury can be classified into mild, moderate, and severe, and repetitive mild traumatic brain injuries are now more commonly referred to as chronic traumatic encephalopathy [75]. While computational models can help estimate injury risk and advance injury prevention [38], correct predictions of the deformation field during impact critically depends on the accurate representation of the nonlinear, loading-mode-, and region-dependent stressstrain relationship of brain tissue. This is especially challenging since impact loading involves extremely high strain rates. It is not surprising that previous studies have shown that simulation results largely depend on both the choice of boundary conditions, for instance concerning the brainskull-interface, and the choice of the constitutive model for brain tissue, as illustrated in Fig. 42 [184]. If only interested in the instant of the impact, which can occur at strain rates as high as 100/s, we can assume that the hysteresis loop closes, as shown in Fig. 32, and use a finite hyperelastic model. Importantly, however, the corresponding parameters should be calibrated using the unconditioned tissue response at equally high rates. Contrary to slow processes such as brain development or neurodegeneration, during traumatic brain injury, the fluid will not drain and will contribute to the tissue stiffness. In accordance with the considerations in Sect. 4.1.2, the Ogden hyperelastic constitutive model has shown best agreement with experimental results [74,168]. Interestingly, numerical studies further indicate that using a linear viscoelastic model may overestimate the deformations within the brain and, consequently, overestimate the injury Reprinted with permission from [178] Fig. 42 Cumulative maximum principal strains during high-rate impact using different material models for brain tissue (top). Radial-circumferential shear strains resulting from in vivo head rotation (bottom). Reprinted with permission from [184] risk if strain-based injury criteria are used [168]. In turn, a study investigating traumatic spinal cord injury in rats found that the compressible Ogden model with a two-term Prony series [110] achieved good agreement with impactor weightdrop experiments.
Taken together, time-dependent effects are critically important, especially on the small time scales of moderate to severe traumatic brain injury. More important than choosing between an elastic, a viscoelastic, or a poro-viscoelastic model is the use of appropriate high strain rate experiments to accurately calibrate the model. In addition to the correct prediction of tissue strains and stresses through numerical simulations, a major research effort in the coming years will be to accurately determine reliable tolerance criteria for cell death and tissue damage to improve injury diagnosis and advance injury prevention. Another major direction in the coming years will be the exploration of the effects of multiple mild repeated injuries to the head with a quest to better understand how tissue damage accumulates over time and how we can best prevent neurodegeneration in chronic traumatic encephalopathy.

Recommendations
The response of brain tissue-more than any other tissue in our body-is hugely sensitive to the length and time scales during loading. This is mainly a result of the ultrasoft, gellike nature of brain tissue, which implies that biochemical effects on very small length and time scales can have tremendous effects on the overall behavior. At first glance, it appears that previous experimental findings are largely contradictory; in reality, these studies tested different "shades" of brain, as illustrated in Fig. 6.
To capture the highly complex mechanical response of brain tissue, several hyper-, visco-, poro-, or poro-viscoelastic constitutive models have been proposed. In this review, we have systematically introduced and discussed models of different level of complexity. Our goal was to show that the selection of an appropriate model, and with it an appropriate set of material parameters, hugely depends on the application of interest: Slow processes such as brain development may be modeled using finite hyperelasticity, moderately fast processes such as brain surgery are predominated by viscoand poro-elastic effects, and extremely fast processes such as brain injury are dominated by the behavior at very high rates.
Concerning material parameters, we have to keep in mind that parameters calibrated from a single experimental loading mode will likely not represent the physiological response of the brain under complex deformation states occurring in vivo. Similarly, parameters for the simulation of slow processes should be calibrated using the conditioned tissue response, while parameters for fast processes should be based on the unconditioned response, ideally probed at equally high rates. Using the appropriate set of material parameters for a particular application is critical because simulated strains outside the calibrated strain range can easily under-or over-estimate the mechanical effect on the brain, and, with it, safety criteria and injury risk. We conclude that due to the inherent bio-chemo-mechanical complexity of the brain, it is critically important to design experiments hand in hand with the modeling approach, and, ideally, even with the clinical application in mind.

Challenges and Perspectives
We conclude by listing our top ten most important steps that we consider critical to push the frontiers in understanding the mechanics of the human brain: • To date, there is not a single constitutive model with a single set of material parameters capable of correctly predicting the response of brain tissue for small and large strains, as well as fast and slow loading rates. Especially for extreme loading situations, associated with the applications in Sect. 5, this drawback should be tackled in the future, through a combination of well-designed experiments and careful modeling. • A challenge in characterizing the constitutive behavior of the brain is its pronounced microstructural heterogeneity. First steps towards capturing local variations in tissue stiffness have been made by using phenomenological models and calibrating regional material parameters. A true step forward would be to design microstructurally motivated constitutive models, which would inherently capture regional heterogeneities at a finer resolution. • Closely related to the previous point, to date, it is not clear, which microstructural components control macroscopic tissue mechanics. While most microstructural investigations have focused on the cellular componentsthe functional elements of brain tissue-from a mechanics point of view, the role of the extracellular matrix should be investigated more thoroughly. Notably, however, load bearing elements in the extracellular matrix of other soft tissues such as collagen do not play a major role in the brain, which explains its incredible compliance compared to other tissues. • Another open question remains, how the in vivo properties of brain tissue compare to the reported ex vivo parameters. Hard tissues have been extensively studied because their mechanical properties change very slowly and only marginally from live to dead. Soft tissues like tendons, ligaments, or arteries, change more drastically but over the past decade, we have learned how to preserve them appropriately post mortem. Ultrasoft tissues like the brain can undergo drastic and rapid changes, both biochemically and biomechanically. To fully understand their in vivo behavior, it will be key to intelligently combine in vivo and ex vivo experiments to develop and calibrate reliable mechanical models for brain tissue in the future. • In addition to the appropriate choice of the constitutive model and the appropriate evaluation of strains, it is essential to accurately capture the boundary conditions during finite element simulations of brain tissue behavior. Especially in the context of high impact loading, during moderate to severe traumatic brain injury, the brainskull interface and internal structures like the falx will significantly affect the deformation profile, and, with it the strains and stresses across the brain. The appropriate choice of boundary conditions needs to be carefully addressed in the future. • Once we understand the constitutive behavior of healthy brain tissue, a next critical step is to predict the risk of axonal injury and tissue damage, and, ultimately, the progression of neurodegenerative disease. It will be essential to perform well-designed experiments to quantify the time course of events during neurodegeneration, and to identify reliable failure criteria under multiaxial loading conditions. Again, these failure criteria will likely depend on the loading mode, as well as the time and length scale of loading. • An important phenomenon that has not been widely considered yet is the fact that brain tissue is not passive but that it contains living cells that actively respond to mechanical stimuli. Brain cells sense and respond to their mechanical environment. Understanding the role of mechanosensing and mechanotransduction is critical in directing migration paths of nervous cells or tumor cells within nervous tissue and in guiding the remodeling of the tissue microstructure. • Understanding the mechanobiology of brain tissue requires a deep biochemical knowledge of the mechanical stimuli that can modulate the response of the brain. Recent studies have highlighted the importance of mechanical stiffness as a guiding principle in neurodevelopment and neurodegeneration. Manipulating the mechanical microenvironment and recreating developmental conditions could hold the key to enhance axonal sprouting, trigger axonal regrowth, and simulate remyelination and provide a new path towards new treatment strategies for neuronal regeneration. • To achieve the long-term goal of using mechanical models for clinical applications, a major challenge is the numerical implementation of sophisticated, highly coupled multiscale and multiphysics problems. When simu-lating whole brain phenomena, the design of efficient and robust algorithms will be critical to make these models applicable and useful for the community at large. • Finally, validating the models, parameters, and simulations is the most critical step in building confidence in the computational modeling of the brain. A validation across the scales, both in space and time, will require concerted efforts of biologists, engineers, and clinical researchers. Large data sets of healthy and diseased brains have already been collected and are being made available to the broad public. Extracting valuable information from these data sets will likely require tools of machine learning to extract valuable features and correlations in space and time. Rather than using machine learning tools as a black box, without any prior information, our mechanical understanding of the brain can define important constraints for these tools to reduce the parameter space and ensure that the results are physically and physiologically meaningful.
With the immense progress in biomedicine, imaging, and computation, we now have the tools to address these frontiers. If we succeed in tackling these challenges, mechanical modeling and computational simulation of the brain will be highly valuable to advance everyday clinical practice by assisting diagnosis and treatment of neurological diseases.