## 2 Materials and Methods
### Anisotropic Stochastic Cell Rescaling
By substituting Equation (5) in Equation (3), we could obtain \(\mathbf{D}\):  
\[D_{ai\beta j}=\frac{\beta_{T}k_{B}T}{3V\tau_{p}}\delta_{ai\beta}\sum_{\gamma}h _{\gamma i}h_{\gamma j}\,. \tag{6}\]  
By direct substitution, and assuming \(U_{\text{el}}=0\), one could obtain the following equation of motion for the lattice vectors:  
\[\text{d}h_{ai}=-\frac{\beta_{T}}{3\tau_{p}}\bigg{[}\sum_{\beta}\Big{(}P_{0} \,\delta_{ai\beta}-P_{\text{int}\text{A}\beta}\Big{)}h_{\beta i}-\frac{k_{B}T} {V}h_{ai}\bigg{]}\,\text{d}t+\sqrt{\frac{2\beta_{T}k_{B}T}{3V\tau_{p}}}\sum_{ \beta}h_{\beta i}\,\text{d}W_{a\beta}\,, \tag{7}\]  
or, in matrix notation:  
\[\text{d}\mathbf{h}=-\frac{\beta_{T}}{3\tau_{p}}\Big{[}\big{(}P_{0}\,\mathbf{ I}-\mathbf{P}_{\text{int}}\big{)}-\frac{k_{B}T}{V}\mathbf{I}\Big{]}\,\mathbf{h}\, \text{d}t+\sqrt{\frac{2\beta\tau k_{B}T}{3V\tau_{p}}}\text{d}\mathbf{W}\, \mathbf{h}\,. \tag{8}\]  
Here, \(\mathbf{P}_{\text{int}}\) is the internal pressure tensor defined as:  
\[P_{\text{int},a\beta}=-\frac{k_{B}T}{V}\sum_{i}h_{\beta i}\frac{\partial(K+U) }{\partial h_{ai}}\,. \tag{9}\]  
As discussed in the next Subsection, the exact interpretation of these partial derivatives depends on how positions and momenta are affected by the changes in the lattice vectors. In presence of external anisotropic stress, Equation (8) had to be amended with the additional contribution associated with the \(U_{\text{el}}\) term, leading to:  
\[\text{d}\mathbf{h}=-\frac{\beta_{T}}{3\tau_{p}}\Big{[}\big{(}P_{0}\mathbf{I}- \mathbf{P}_{\text{int}}\big{)}-\frac{k_{B}T}{V}\mathbf{I}+\frac{\mathbf{h} \mathbf{\Sigma}\mathbf{h}^{T}}{V}\Big{]}\,\mathbf{h}\,\text{d}t+\sqrt{\frac{2 \beta_{T}k_{B}T}{3V\tau_{p}}}\text{d}\mathbf{W}\,\mathbf{h}\,. \tag{10}\]  
Equation (8) (or (10)) had to then be coupled to the propagation of the Hamilton equations for positions and velocities, at constant cell matrix, and to a thermostat so as to control temperature in addition to pressure. Any thermostat could be used, but one should take into account that, if the potential energy \(U\) is invariant with respect to a rigid translation of the coordinates and a global thermostat is used, the velocity of the center of mass would be constrained. To make the results from a global thermostat equivalent to those obtained with a local thermostat, one might need to add the center of mass contribution to the internal pressure defined in Equation (9) [12; 18; 23]. In the current work, we only tested our barostat in combination with a global thermostat [24] and decided not to include this correction so as to obtain results equivalent to those obtained with the commonly used MTTK algorithm [8; 9]. Equation (8) can be shown to be invariant with respect to linear combinations of the lattice vectors, so that the resulting dynamics do not depend on the arbitrary choice of lattice vectors associated with a given Bravais lattice. In addition, it can be proved that, if the scaling is constrained to be isotropic or semi-isotropic, the scheme reduces to the isotropic or semi-isotropic versions of the SCR scheme [18], respectively. The proofs of these statements are sketched in Appendix A, while the full derivations were reported in Ref. [25]. Furthermore, it is possible to show that Equation (10) could be obtained from a Parrinello-Rahman barostat augmented with friction and noise by taking the limit of infinite friction tensor \(\gamma\) and null inertia \(W\) at fixed \(\gamma W\). As discussed in detail in Ref. [25], this derivation brings to Equation (10) all but a correction term that is negligible in the thermodynamic limit, originating from the absence of the factor \((\det\mathbf{h})^{-2}\) in the definition of the \(N\mathbf{S}T\) ensemble used in the Parrinello-Rahman method. Even though the Parrinello-Rahman method is not invariant with respect to linear combinations of the lattice vectors [7], the resulting SCR scheme is invariant thanks to the form adopted for the tensor \(\gamma W\)[25].