Modelling Three-Dimensional Microstructure Evolution Influenced by Concurrent Structural Mechanical Mechanisms

The interdependence between structural mechanics and microstructure solidification is an inherently three-dimensional phenomenon, where the complex physical processes and mechanical interactions can lead to dendrites growing at orientations influenced by twisting and out of plane bending. These effects can have a significant impact on the formation of defects and the overall macroscopic material properties of the structure. However, all attempts to numerically model this process so far have been limited to two-dimensional representations of the problem, which necessitates ignoring any potential behaviour that may arise from these more complex deformation events. For this reason, the two-dimensional numerical methods presented in previous papers, which couple a Finite Volume Structural Mechanics Solver to a Cellular Automata solidification solver, have been expanded so that problems may now be simulated in three dimensions. Results are presented which do not aim to predict any specific mechanism but rather highlight the new capabilities of this improved three-dimensional modelling framework.


INTRODUCTION
The impact of structural mechanics upon the microstructural development of solidifying metal alloys is currently an under-explored field of research, particularly given the changes these effects have been observed to have on dendrite development and composition, dictating the macroscopic material properties of the ultimate solidified structure. The interdependent nature of concurrent structural mechanics and microstructure solidification can lead to complex interactions, where any forces which may cause mechanical deformation will necessarily alter subsequent solidification. In extreme cases large deformations can lead to dendrite fragmentation 1 and act as a contributing factor in the formation of structural inhomogeneities. 2 In more extreme cases other casting defects such as stray grains 3,4 and slivers 5-7 may arise, which necessitate parts being discarded as unfit for purpose.
A key factor in the interaction between structural mechanics and microstructure solidification is the crystallographic orientation of the dendrite, which dictates how it will preferentially grow. As a dendrite deforms the local orientation throughout the dendrite will correspondingly change, influencing the continuous development of the microstructure. Changes to dendrite orientation can be gradual, with successive small misorientations accumulating along a dendrite arm until a significant overall change to the orientation is achieved. 8,9 However there are also cases of large localized orientation changes, where dendrite growth will continue at the new orientation with no further observable orientation changes. 10,11 Due to the complex nature of these interactions, it can simplify matters for both modelling and practical experiments to treat systems as being functionally two dimensional. For practical experiments, due to limitations in obtaining in situ visuals of the microstructure as it develops, thin samples can be utilized where one of the dimensions takes a This article was updated to replace https://doi.org/10.1007/s11837-022-05232-0 Ó 2022 The Author(s) negligible size to allow for a complete observation of the experiment. 10,12,13 In cases where larger samples are used, analysis is generally relegated to post processing the final solidified structure to examine any potential regions of interest, with the available methods generally necessitating that any appraisal of the microstructure is performed on two-dimensional images describing slices of the microstructure, 1,3,5,11 which can complicate both identifying and visualizing more complex mechanical behaviours such as torsion and out-of-plane bending. 10,14 Numerically modelling concurrent structural mechanics and microstructure solidification could provide a greater understanding of how these interactions change the transient development, not being bound by the same constraints as practical experiments when it comes to obtaining in situ visuals. Furthermore, with a fully simulated dendritic structure there is potential for the various material properties of the fully solidified macrostructure to be predicted. However, numerical modelling of these effects remains at a nascent stage with structural mechanics generally only being applied to post process macroscopic structures, [15][16][17][18] overlooking any impact structural mechanics will have on further solidification.
Numerical modelling has been undertaken which considers the impact of structural mechanics on the microstructure, in the context of both post processing 19,20 and fully coupled simulations. 21 However, all existing modelling considers the dendrites only in two dimensions where a single columnar dendrite is examined and assumed to behave in a manner analogous to a cantilever beam. This two-dimensional approach was also adopted in previous work published by the authors, with Soar et al. 22 investigating the influence to the growth behaviour from applying forces which complement and oppose the growth direction. This was extended in Soar et al. 23 to model larger dendritic systems compared to experimental results obtained using a Ga-In alloy, where forces perpendicular to the growth orientation of the primary dendrite arm were considered as a source of dendrite misorientation which may influence grain competition.
However, while a two-dimensional model may act as a fair representation for thin sample cases for many scenarios, even within this framework the process is intrinsically three-dimensional in nature, where the microstructure may develop in ways twodimensional simulations could not capture. Consequently, modelling behaviour such as dendrites bending out of plane, twisting, or overlapping in the sample is rendered impossible by remaining in two dimensions.
The modelling presented by Soar et al. 23 has been expanded so that the microstructural evolution is modelled in three dimensions, with the intimately coupled structural mechanical behaviour and locally defined crystallographic orientation likewise expanded to function in three dimensions. This expansion of coupled structural mechanics and microstructure solidification to 3D allows the examination of geometric changes which cannot be represented using any other currently existing numerical models which are limited by working in two dimensions. The modelling is currently limited to examinations of columnar dendrites growing under directional solidification conditions typical for casting. While this approach can model an equiaxed dendrite bound by contact with other grains, if the dendrite is free floating it will be unbound from a structural mechanical perspective, requiring a further solid body model to describe their movement through the domain. Results have been generated which examine idealized cases which highlight the capabilities of this unique modelling framework.

MATHEMATICAL FORMULATION
The system being solved using the numerical method presented in Soar et al. 23 has been expanded to represent a three-dimensional system of concurrent microstructure solidification and structural mechanics.
The microstructure evolution is unchanged from prior work by being represented in three dimensions, with the governing equations for phase change and partitioning being given by where C l and C s are respectively the liquid and solid concentration, k ¼ C l =C s is the partition coefficient, 0 / s 1 is the solid fraction, D e is the mass diffusivity, and t is time. Diffusive transport is dictated by the relationship with the equivalent concentration defined as The equilibrium temperature T i is given by where m l is the liquid slope, C 0 is the initial concentration, and T 0 is the liquidus temperature at C 0 . The structural mechanical behaviour is represented by a linear elastic material model, with the corresponding assumptions of 'small' deformations, no plastic deformation and a linear relationship between the stress and strain. In modelling the third dimension, the assumption of plane strain used in earlier work can now be discarded. In three dimensions, to obtain the stationary solution for the displacements of a structure experiencing external forces, the following partial differential equations need to be simultaneously solved: where u, v and w are the displacements in the x, y and z direction, F x , F y and F z are the corresponding body forces, and l and k are the Lamé constants: where E and g are respectively the Young's modulus and Poisson's ratio of the material. Variable material properties are accounted for by defining the Young's modulus as a parameter which varies locally according to the following linear relationship with the solid fraction and the maximum Young's modulus ðE M Þ of the fully solidified material: In the cases presented, a fixed body force, F ¼ ðq s À q l Þg, based the density difference between the solid (q s Þ and liquid (q l Þ, is applied to the solidifying structure. F can be defined as acting in any combination of the x, y and z directions as the problem necessitates and g acts as a gravitational constant which takes multiple sizes along with both positive and negative signs depending on the case being examined, as will be explicitly stated in the problem descriptions.

NUMERICAL METHOD
The elastic deformation of dendrite arms is computed using a locally defined variable, Young's modulus, as a parameter, changing crystallographic orientation of the growing dendrites. In the absence of experimental data, this was proposed as an approximate model of the observed complex elastovisco-plastic phenomena occurring at the microscale during alloy solidification under external forces or geometrical constraints.
The microstructure evolution is modelled using a finite-difference decentred octahedral Cellular Automata (CA) method. This method is based on the lMatIC code [24][25][26][27] to resolve the evolution of the alloy microstructure, which has been further developed by the authors to be able to solve complex multi-physics phenomena in parallel within the modelling framework TESA 28,29 where it has been extensively verified and validated. This code was already capable of solving both two-and threedimensional problems prior to the development of the Structural Mechanics Solver (SMS), so no further development relating to the microstructure solidification process was required.
The structural mechanics is solved using a bespoke Finite Volume Structural Mechanics Solver which has been developed to solve the linear elasticity equations in displacement formulation, allowing for easy coupling with the CA method used by TESA. While applying a finite volume approach to structural mechanics is somewhat unusual, a significant body of work exists which has rigorously verified and validated this approach, finding it to be comparable with solutions obtained using the finite element method. [30][31][32] Other than requiring a new discretisation to solve these updated equations in three dimensions, the implementation of the SMS remains the same as was described in Soar et al., 23 remaining capable of obtaining the displacements for structures with variable material properties. This is treated as a quasi-stationary process where only new displacements arising since the prior call of the solver are calculated, requiring the body forces to act transiently such that they only apply the change in force since the prior timestep, rather than the force in total.
The same process of coupling structural mechanics to microstructure solidification is employed by means of locally changing the crystallographic orientation of a dendrite based upon the deformation it has experienced. The primary way this process is changed by working in three dimensions is that while the dendrite orientation could be described by a single angle in two dimensions, in three-dimensional space it is now required to calculate three extrinsic rotations around the axes to describe all possible orientations a dendrite may take. In practice, obtaining these rotations leads in effect to separate 2D problems where the rotation about each axis can be calculated by using the obtained displacements to form an arc as described in Soar et al. 23 and visualised in Fig. 1. This approach allows the three extrinsic rotations around the x, y and z to be calculated as follows: Dh y ¼ 1 2Dx u h À u l À w e þ w w ð Þ % 1 2 @u @z À @w @u ð11Þ where the subscripts e, w, n, s, h and l indicate that the displacements being used are spatially located respectively at the East, West, North, South, High and Low faces of the volume whose orientation is being calculated. Once obtained, these orientation changes can be combined with any existing orientation information at that point to obtain new values of h x ; h y and h z which can be used to generate a rotation matrix R ¼ R z h z ð ÞR y h y À Á R x h x ð Þ for each volume which the CA method interprets to describe the local orientations dictating how dendrites preferentially grow in a three-dimensional space.
In summary, for a single transient solution step using these methods, the CA method will provide the evolution of the microstructure that will dictate changes to the material properties and body force being applied. This information is then passed to the SMS which solves the linear elasticity equations in three dimensions to obtain displacement values which are used to alter the local crystallographic orientation of each volume comprising the structure. The CA method will then be called again, considering the updated three-dimensional orientation of the structure as the solution progresses.

RESULTS
The results presented have been chosen to highlight the capabilities of the updated numerical model, rather than necessarily predict any specific mechanism. The first case chosen is a thin sample case analogous to a two-dimensional problem in expected behaviour, because one of the dimensions of the sample goes into the plane of view is significantly smaller than the other dimensions of the sample domain. For the purposes of modelling thin sample experiments using the now 3D SMS, a sample width of 200lm was used, meaning that the behaviour observed should be comparable to those seen in a two-dimensional case being modelled using plane stress. The relevant material properties for the majority of the cases presented are given in Table I, being analogous to those of a semi-solid Ga-25 wt%. In alloy, chosen for easier comparison with the prior work published by Soar et al. 23 (with a lower value taken for the maximum Young's modulus to highlight the bending behaviour in the smaller simulation domains) which modelled a system using this same alloy to allow for comparison with experimental results. All cases presented use a cell size of Dx ¼ 10lm, a timestep Dt ¼ 5ms, a thermal gradient of rT ¼ 1K=mm and a cooling rate Q ¼ 0:01K=s. This leads to characteristic solidification rates of 10lm=s, which correspond to typical conditions for casting. External body forces have been applied in all cases to cause observable deformation, which highlights the capabilities of the SMS, which may come from a wide range of sources depending on the manufacturing process. This could include pressures from fluid flow, thermal contraction/shrinkage, electromagnetic forces or density variations. For the cases being presented, forces have been selected that are simple to define, being based on a constant gravitational acceleration acting on the inherent density variations caused by partitioning, or a fixed rotational force to demonstrate the three-dimensional behaviours of the  coupled system. In principle, the model can account for any external body force, but the application of many of these would require additional processes to be modelled to obtain the requisite source values. The first case presented in Fig. 2 takes a solidifying single columnar dendrite while the structure is subjected to a transient body force based on the density difference of the growing structure and surrounding solute. The setup for this case takes a single seed placed on the west wall of a 400 9 160 9 20-cell domain with an initial orientation of 5°, taking just under 10 h to run on a 16-core cluster node. In Fig. 2(a-c) the transient evolution of the dendrite can be observed as it grows across the domain, where in (a) only a slight deformation is beginning to be visible; by (b) the orientation has changed to the degree that the starting orientation of 5°has almost been entirely counteracted while by (c) the dendrite now clearly has an orientation heading in the opposite direction to the angle the dendrite was initially seeded with. These orientation changes are highlighted further in (d), which shows the change in orientation growing along the length of the dendrite, with the largest orientation change being observed at the tip. Some small positive changes in orientation can be observed in the secondary arms growing from the wall, being caused by the elastic nature of the material which is causing this region to bend upwards to respond to the deformation in the main arm. Part (e) shows the von Mises stress within the dendrite, highlighting that the highest regions of stress are generally where secondary arms attach to the main arm or where they interact with each other. Finally, (f) shows the state of a dendrite grown under the same conditions without a force being applied for the sake of comparison. The results in (c) and (f) have been presented as a 3D structure rather than slices through it to demonstrate the true three-dimensional morphology of the structures being modelled.
This test was then expanded to model a second case of a larger thin sample, comparable to the 2D simulation presented in Soar et al. 23 This took a 6400 9 1600 9 20-cell domain where two dendrites were seeded on the west wall with an initial 20°o rientation, where the North and South boundaries are periodic such that this system could be considered as representative of a section of a larger sample. This case was run on the Greenwich University HPC, taking approximately 16 h to complete utilising 400 cores. Contrary to the other results presented in this paper, this case was modelled subject to the material properties outlined in the original 2D simulation 23 rather than using the values given in Table I. Consequently, the results shown in Fig. 3 were obtained by modelling conditions taking E M = 30 MPa and the body force being equivalent to three times the terrestrial gravity. These results demonstrate similar behaviour overall to that observed in the earlier 2D simulation, but this scenario has also captured a new behaviour not observed in earlier modelling. This can be observed in (a) where some of the primary arms are bending away from each other due to secondary arm interaction, which if continued could lead to a new primary arm forming in the space growing between them and the competition of the converging arm tips leading to other dendrites becoming trapped. Part (b) shows the same behaviour with accumulated deformation as was observed in the 2D simulation is in evidence here, having a maximum in the centre of the dendrite with the values fading towards the root and tip of the dendrites. Only the v deformations in y have been presented, as the dominant deformation values driving orientation change in this simulation. The orientations around the z axis shown in (c) serve to better highlight the beginning divergence of the dendrite arms. For this case orientation changes around the x and y axes have also been calculated, but the values were found to be insignificant in affecting the overall solution. The von Mises stress in (d) shows stress accumulating throughout the dendritic structure except for the tips, which are still relatively free to move, with there being a high region of stress that quickly diminishes at the point where the arms start to noticeably diverge. It can also be observed that the deformations and corresponding orientation changes are significantly lower than those found in the 2D simulations, despite using the same modelling conditions. The reason for this inconsistency is that two-and threedimensional solidification are inherently different, with an extra degree of freedom for partitioning and mass transport in 3D causing dendrites to essentially become thicker, as can be seen in analytic solutions. 33 These thicker dendrites will be correspondingly harder to bend, making direct comparisons between 2D and 3D an involved process as it is not just boundary condition changes that need to be considered, but the entire system as phase change and transport play a role. Consequently, the main conclusion to be drawn from this comparison is that within a thin sample context this approach can obtain fully 3D results for domain sizes comparable to those used in practical experiments which were previously modelled only in 2D, where this added dimension can significantly alter the observed behaviour.
Better examples to demonstrate behaviour in three dimensions are cases where the loading causes a dendrite to bend out of plane. While this can certainly occur in a thin sample case, it is more common and obvious in fully 3D cases where all dimensions of the domain have a significant size. Consequently, for the next case a 120 9 120 9 480cell domain with periodic conditions on the East, West, North and South boundaries was modelled, with the simulation taking approximately 11 h using a four-core machine. Here a single seed with no pre-defined orientation was placed on the Low boundary and allowed to develop under negative body forces acting in the x and y directions, both taking a g value five times larger than terrestrial gravity. This caused the columnar dendrite to bend into one of the corners of the domain as can be observed in Fig. 4. A view of the dendrite has been provided from multiple angles while displaying the total local orientation of the dendrite around the x and y axes. The local orientation changes around the z axis were also calculated during the simulation, but these values are insignificant in magnitude compared to the other orientations.
With the capability of the SMS to capture more complex three-dimensional bending demonstrated by the previous case, a natural extension would be to examine how these structural effects can change the development of fully three-dimensional samples featuring multiple dendrites. This has been realized in the case presented in Fig. 5, where 50 dendrites were seeded at a random location with a random starting orientation on the low face of a 120 9 120 9 480-cell domain with periodic conditions on the East, West, North and South boundaries, with the simulation taking approximately 8 h using a fourcore machine. Thus, this case can be considered representative of a sub-sample of a larger threedimensional domain. This system was then modelled under conditions with no force and under a constant negative body force applied in the z direction taking a g value equivalent to five times terrestrial gravity. This highlights the fact that  within the complex framework of competing dendrites in three dimensions, the SMS can correctly identify and obtain displacements for the individual structures, while also keeping track of the distinct local orientations occurring within these dendrites. While the same two dendrites ultimately outcompete the others to become the primary arms in both cases, clear differences in the dendrite orientation can be observed throughout the sample. For larger samples or different initial conditions based on the random seed used, it would not be unreasonable to expect the impact from these structural mechanical changes to influence which dendrites successfully manage to outcompete their neighbours.
The final case being examined in this paper studies the application of a rotational force to a growing dendrite, with twisting dendrites being another experimentally observed behaviour which can now be represented using the SMS. For this, a 120 9 120 9 480-cell domain with periodic conditions on the East, West, North and South boundaries was again used, with the simulation taking approximately 10 h using a four-core machine. A single seed with no pre-defined orientation was placed on the Low boundary so that it can grow through the rotational force being applied to the domain. This rotational force was represented by taking body forces acting in the x and y directions, which have a value of zero at the centre of the domain, but which increase radially in magnitude such that the g reaches positive or negative values of five times terrestrial gravity at the domain boundaries. This imparts a counterclockwise force upon the structure, causing it to develop into a clear spiral as observed in Fig. 6. In 6(a-c), this shows the time evolution of the structure, which demonstrates that not only do the secondary arms twist such that they preferentially grow counterclockwise, but the central arm itself twists as it grows, shifting the location along the arm where secondary arms begin their growth. Part (d) shows how this rotation around the z axis varies along the length of the dendrite such that the tip has rotated a full 360°d uring the simulation. Notwithstanding the linear elasticity assumption, this demonstrates that even under quite extreme loading conditions, the SMS can track orientation change imparted by torque, allowing for the modelling of scenarios where structures twist as they grow.

CONCLUSION
The two-dimensional structural mechanics code presented in an earlier work has been extended to solve the linear elasticity equations for a threedimensional evolving structure. This necessitated a corresponding expansion of the local orientation calculations, such that rotations around all axes could be obtained from the deformations, which could then be used to describe the orientation changes in three-dimensional space. The updated structural mechanics code was coupled to a CA method solidification code using deformations to update crystallographic orientations. This SMS-CA code was then applied to a selection of cases that highlights the unique capabilities of this modelling framework. Initially, thin sample cases analogous to the two-dimensional solutions previously examined were simulated to verify the code was still capable of resolving this style of domain which is commonly used in practical experiments where three-dimensional behaviour is generally limited. The remaining cases examined behaviours which can only exist in three-dimensions, such as out-of-plane bending and twisting. The cases presented in this paper are idealized to demonstrate the new capabilities of the three-dimensional code which cannot be accounted for by using any other currently existing numerical models. This serves to highlight the potential of this framework to capture complex behaviour observed in experiments which are often associated with defects.

FUTURE WORK
With the three-dimensional implementation complete, the next stage of development would be the integration of the other physics, for example allowing structural mechanics to be solved concurrently with fluid flow. This would alter the solute distribution to change the microstructure development as well as potentially directly impart additional forces upon the structure. Other avenues of further development would be the implementation of a material model capable of resolving non-linear structural mechanics and thermal effects. These developments would allow for the numerical modelling of more realistic cases which can better capture the behaviour observed in practical experiments.