diff --git a/cookbooks/CPO_induced_anisotropic_viscosity/doc/CPO_induced_anisotropic_viscosity.md b/cookbooks/CPO_induced_anisotropic_viscosity/doc/CPO_induced_anisotropic_viscosity.md index 96738311d7c..515f55c4d8d 100644 --- a/cookbooks/CPO_induced_anisotropic_viscosity/doc/CPO_induced_anisotropic_viscosity.md +++ b/cookbooks/CPO_induced_anisotropic_viscosity/doc/CPO_induced_anisotropic_viscosity.md @@ -16,9 +16,9 @@ This cookbook explains how to use the CPO-induced anisotropic viscosity material ## Introduction -Individual crystals of the mineral olivine reorganize their orientations into crystal-preferred orientations (CPO) under deformation. The viscous properties of olivine crystals are direction-dependent (anisotropic), which suggests that the effective viscosity for olivine rocks/aggregates is different when deformations occur in different directions relative to the CPO. This cookbook model computes an anisotropic viscosity based on the CPO evolution predicted by D-Rex ({cite}`fraters_billen_2021_cpo`; {cite}`kaminski2004`) and includes this information in the subsequent modeling process. +Individual crystals of the mineral olivine reorganize their orientations into crystal-preferred orientations (CPO) under deformation. The viscous properties of olivine crystals are direction-dependent (anisotropic), which suggests that the effective viscosity for olivine rocks/aggregates is different when deformations occur in different directions relative to the CPO. This cookbook model computes an anisotropic viscosity based on the CPO evolution predicted by D-Rex ({cite}`fraters_billen_2021_cpo`; {cite}`kaminski2004`) and includes this information in the subsequent modeling process based on the formulation of Király at al., (in rev.). -Our constitutive equation for the relationship between the strain rate and stress using the anisotropic viscosity tensor is adapted from {cite:t}`signorelli:etal:2021`: +The constitutive equation for the relationship between the strain rate and stress using the anisotropic viscosity tensor is adapted from {cite:t}`signorelli:etal:2021`: ```{math} :label: eqn:anisotropic_general_stress @@ -36,7 +36,7 @@ $\gamma_0=1.1\times 10^{5}$ is the isotropic fluidity, $Q=530$ $kJ/mol$ is the a ```{math} :label: eqn:equivalent_yield_stress -J(\sigma_{ij})=(F(\sigma_{11} - \sigma_{22})^2+G(\sigma_{22} - \sigma_{33})^2+H(\sigma_{33} - \sigma_{11})^2+2L\sigma_{23}^2+2M\sigma_{13}^2+2N\sigma_{12}^2)^{1/2} +J(\sigma_{ij})=(F(\sigma_{22} - \sigma_{33})^2+G(\sigma_{33} - \sigma_{11})^2+H(\sigma_{11} - \sigma_{22})^2+2L\sigma_{23}^2+2M\sigma_{13}^2+2N\sigma_{12}^2)^{1/2} ``` and $A_{ij}$ is the anisotropic tensor of fluidity in Kelvin notation: @@ -45,24 +45,24 @@ and $A_{ij}$ is the anisotropic tensor of fluidity in Kelvin notation: :label: eqn:anisotropic_fluidity A_{ij}=\frac{2}{3} \left[ \begin{matrix} -F+H & -F & -H & 0 & 0 & 0 \\ --F & G+F & -G & 0 & 0 & 0 \\ --H & -G & H+G & 0 & 0 & 0 \\ +G+H & -H & -G & 0 & 0 & 0 \\ +-H & H+F & -F & 0 & 0 & 0 \\ +-G & -F & F+G & 0 & 0 & 0 \\ 0 & 0 & 0 & L & 0 & 0 \\ 0 & 0 & 0 & 0 & M & 0 \\ 0 & 0 & 0 & 0 & 0 & N \end{matrix} \right] ``` -$J(\sigma_{ij})$ and $A_{ij}$ are computed using Hill coefficients $H, J, K, L, M,$ and $N$ {cite}`hill:1948`, which describe the anisotropic viscous properties of an olivine aggregate and depend on its CPO. We determine the mean CPO orientation from the eigenvectors associated with the largest eigenvalues of the second-order orientation tensor (or covariance matrix) for all three symmetry axes. The corresponding eigenvalues quantify the dispersion of orientations around these mean orientation {cite}`bingham:1974`. The relationship between the 9 eigenvalues (3 for each axis) and Hill coefficients is derived using regression analysis on a texture database constructed with olivine textures from laboratory experiments, shear box models, and subduction models (Kiraly et al., in rev.). The 9 coefficients and 1 constant for each of the Hill coefficients are given as input in the parameter file. The default values and the equation to compute the Hill coefficients from the eigenvalues (e.g. $a_1, a_2, a_3$ are the eigenvalues of the orientation tensor for a-axis, where $a_1$ is the largest eigen value) are shown below: +$J(\sigma_{ij})$ and $A_{ij}$ are computed using Hill coefficients $H, J, K, L, M,$ and $N$ {cite}`hill:1948`, which describe the anisotropic viscous properties of an olivine aggregate and depend on its CPO. We determine the mean CPO orientation from the eigenvectors associated with the largest eigenvalues of the second-order orientation tensor (or covariance matrix) for all three symmetry axes. The corresponding eigenvalues quantify the dispersion of orientations around these mean orientation {cite}`bingham:1974`. The relationship between the 9 eigenvalues (3 for each axis) and Hill coefficients is derived using regression analysis on a texture database constructed with olivine textures from laboratory experiments, shear box models, and subduction models (Kiraly et al., in rev., the data used for this paper is published through Zenodo {cite}`kiraly:etal:2026`). The 9 coefficients and 1 constant for each of the Hill coefficients are given as input in the parameter file. The default values and the equation to compute the Hill coefficients from the eigenvalues (e.g. $a_1, a_2, a_3$ are the eigenvalues of the orientation tensor for a-axis, where $a_1$ is the largest eigen value) are shown below: ```{math} :label: eqn:hill_coefficients -F = 1.039 a_1^2 - 0.767 a_2 - \frac{0.003}{a_3} + 0.197 b_1^2 + 0.413 b_2 + \frac{0.015}{b_3} - 0.936 c_1^2 - 2.393 c_2 + \frac{0.052}{c_3} + 1.08 \\ -G = -2.836 a_1^2 - 1.632 a_2 - \frac{0.001}{a_3} + 0.267 b_1^2 - 0.993 b_2 + \frac{0.003}{b_3} +1.969 c_1^2 + 2.314 c_2 - \frac{0.019}{c_3} + 0.69 \\ -H = 1.669 a_1^2 + 0.58 a_2 + \frac{0.003}{a_3} + 0.702 b_1^2 + 0.251 b_2 + \frac{0.000}{b_3} - 2.003 c_1^2 - 2.570 c_2 + \frac{0.071}{c_3} + 0.75 \\ -L = -0.325 a_1^2 + 0.728 a_2 + \frac{0.000}{a_3} - 0.665 b_1^2 + 0.515 b_2 + \frac{0.003}{b_3} - 1.027 c_1^2 - 1.263 c_2 + \frac{0.009}{c_3} + 1.60 \\ -M = 1.643 a_1^2 - 0.878 a_2 + \frac{0.005}{a_3} + 2.489 b_1^2 + 0.816 b_2 - \frac{0.011}{b_3} - 2.494 c_1^2 - 0.511 c_2 + \frac{0.009}{c_3} + 0.89 \\ -N = 0.812 a_1^2 - 0.157 a_2 + \frac{0.002}{a_3} - 1.649 b_1^2 + 0.194 b_2 - \frac{0.01}{b_3} + 1.68 c_1^2 - 0.104 c_2 + \frac{0.02}{c_3} + 1.21 +F = 0.592 a_1^2 - 0.832 a_1 - 0.001 a_2 - \frac{0.000}{a_3} + 0.380 b_1^2 - 0.533 b_1 + 0.468 b_2 - \frac{0.001}{b_3} - 1.249 c_1^2 + 1.075 c_1 - 0.168 c_2 + \frac{0.003}{c_3} + 0.52 \\ +G = -1.695 a_1^2 + 1.336 a_1 - 0.184 a_2 + \frac{0.000}{a_3} + 0.750 b_1^2 + 0.691 b_1 + 0.377 b_2 - \frac{0.002}{b_3} - 0.670 c_1^2 - 0.552 c_1 - 0.428 c_2 + \frac{0.003}{c_3} + 0.26 \\ +H = -1.140 a_1^2 + 1.353 a_1 + 0.751 a_2 - \frac{0.002}{a_3} - 0.256 b_1^2 - 1.006 b_1 - 0.116 b_2 + \frac{0.003}{b_3} + 0.648 c_1^2 - 0.031 c_1 - 0.080 c_2 + \frac{0.006}{c_3} + 0.75 \\ +L = -3.511 a_1^2 + 2.686 a_1 + 0.360 a_2 - \frac{0.001}{a_3} + 3.948 b_1^2 - 3.816 b_1 - 0.779 b_2 + \frac{0.004}{b_3} + 4.122 c_1^2 - 2.483 c_1 - 1.320 c_2 + \frac{0.002}{c_3} + 2.00 \\ +M = 4.537 a_1^2 - 3.228 a_1 + 0.276 a_2 + \frac{0.007}{a_3} - 7.447 b_1^2 + 5.764 b_1 - 1.403 b_2 - \frac{0.032}{b_3} + 2.968 c_1^2 - 3.435 c_1 - 2.266 c_2 + \frac{0.122}{c_3} + 2.44 \\ +N = 7.873 a_1^2 - 7.934 a_1 - 2.588 a_2 + \frac{0.030}{a_3} + 7.606 b_1^2 - 5.469 b_1 - 0.348 b_2 + \frac{0.064}{b_3} - 1.788 c_1^2 + 2.255 c_1 + 3.023 c_2 - \frac{0.103}{c_3} + 3.70 ``` In this material model plugin, strain rate, density, temperature, and other parameters are taken as input to compute the anisotropic viscosity, which is passed into the Stokes system to compute the stress. As a result, we adapt {math:numref}`eqn:anisotropic_general_stress` to be: diff --git a/cookbooks/CPO_induced_anisotropic_viscosity/plugin/CMakeLists.txt b/cookbooks/CPO_induced_anisotropic_viscosity/plugin/CMakeLists.txt index 38bce63c4d5..6dd63c7b2c1 100644 --- a/cookbooks/CPO_induced_anisotropic_viscosity/plugin/CMakeLists.txt +++ b/cookbooks/CPO_induced_anisotropic_viscosity/plugin/CMakeLists.txt @@ -16,24 +16,25 @@ # along with ASPECT; see the file LICENSE. If not see # . -CMAKE_MINIMUM_REQUIRED(VERSION 3.13.4) +cmake_minimum_required(VERSION 3.13.4) -FIND_PACKAGE(Aspect 2.4.0 QUIET HINTS ${Aspect_DIR} ../ ../../ $ENV{ASPECT_DIR}) +find_package(Aspect 2.4.0 QUIET HINTS ${Aspect_DIR} ../ ../../ $ENV{ASPECT_DIR}) -IF (NOT Aspect_FOUND) - MESSAGE(FATAL_ERROR "\n" - "Could not find a valid ASPECT build/installation directory. " - "Please specify the directory where you are building ASPECT by passing\n" - " -D Aspect_DIR=\n" - "to cmake or by setting the environment variable ASPECT_DIR in your shell " - "before calling cmake. See the section 'How to write a plugin' in the " - "manual for more information.") -ENDIF () +if (NOT Aspect_FOUND) + message(FATAL_ERROR "\n" + "Could not find a valid ASPECT build/installation directory. " + "Please specify the directory where you are building ASPECT by passing\n" + " -D Aspect_DIR=\n" + "to cmake or by setting the environment variable ASPECT_DIR in your shell " + "before calling cmake. See the section 'How to write a plugin' in the " + "manual for more information.") +endif () DEAL_II_INITIALIZE_CACHED_VARIABLES() -SET(TARGET "CPO_induced_anisotropic_viscosity") -PROJECT(${TARGET}) -ADD_LIBRARY(${TARGET} SHARED cpo_induced_anisotropic_viscosity.h cpo_induced_anisotropic_viscosity.cc anisotropic_stress.cc anisotropic_stress.h) +set(TARGET "CPO_induced_anisotropic_viscosity") +project(${TARGET}) + +add_library(${TARGET} SHARED cpo_induced_anisotropic_viscosity.h cpo_induced_anisotropic_viscosity.cc anisotropic_stress.cc anisotropic_stress.h) ASPECT_SETUP_PLUGIN(${TARGET}) diff --git a/cookbooks/CPO_induced_anisotropic_viscosity/plugin/cpo_induced_anisotropic_viscosity.cc b/cookbooks/CPO_induced_anisotropic_viscosity/plugin/cpo_induced_anisotropic_viscosity.cc index 979eadff8a1..38b3bc7c7ad 100644 --- a/cookbooks/CPO_induced_anisotropic_viscosity/plugin/cpo_induced_anisotropic_viscosity.cc +++ b/cookbooks/CPO_induced_anisotropic_viscosity/plugin/cpo_induced_anisotropic_viscosity.cc @@ -191,7 +191,7 @@ namespace aspect { // Create constant value to use for AV const double A_o = 1.1e5*std::exp(-530000/(8.314*in.temperature[q])); - const double n = 3.5; + const double n = 3.5; //n=3 for test against VPSC, n=3.5 for D-Rex in ASPECT // The values of A_o and 0.73 were picked so that Gamma = 3.5322e-15[1/(s*Pa^n)] if T=1600K and d=1000 microns const double Gamma = (A_o/(std::pow(grain_size,0.73))); @@ -212,16 +212,16 @@ namespace aspect const double theta = composition[cpo_bingham_avg_b[0]]; const double phi2 = composition[cpo_bingham_avg_c[0]]; - const Tensor<2,3> R = transpose(euler_angles_to_rotation_matrix(phi1, theta, phi2)); + const Tensor<2,3> R = euler_angles_to_rotation_matrix(phi1, theta, phi2); // Compute Hill Parameters FGHLMN from the eigenvalues of a,b,c axis - // CPO2Hill v5 model: - const double F = Utilities::fixed_power<2>(eigvalue_a1)*CnI_F[0] + eigvalue_a2*CnI_F[1] + (1/eigvalue_a3)*CnI_F[2] + Utilities::fixed_power<2>(eigvalue_b1)*CnI_F[3] + eigvalue_b2*CnI_F[4] + (1/eigvalue_b3)*CnI_F[5] + Utilities::fixed_power<2>(eigvalue_c1)*CnI_F[6] + eigvalue_c2*CnI_F[7] + (1/eigvalue_c3)*CnI_F[8] + CnI_F[9]; - const double G = Utilities::fixed_power<2>(eigvalue_a1)*CnI_G[0] + eigvalue_a2*CnI_G[1] + (1/eigvalue_a3)*CnI_G[2] + Utilities::fixed_power<2>(eigvalue_b1)*CnI_G[3] + eigvalue_b2*CnI_G[4] + (1/eigvalue_b3)*CnI_G[5] + Utilities::fixed_power<2>(eigvalue_c1)*CnI_G[6] + eigvalue_c2*CnI_G[7] + (1/eigvalue_c3)*CnI_G[8] + CnI_G[9]; - const double H = Utilities::fixed_power<2>(eigvalue_a1)*CnI_H[0] + eigvalue_a2*CnI_H[1] + (1/eigvalue_a3)*CnI_H[2] + Utilities::fixed_power<2>(eigvalue_b1)*CnI_H[3] + eigvalue_b2*CnI_H[4] + (1/eigvalue_b3)*CnI_H[5] + Utilities::fixed_power<2>(eigvalue_c1)*CnI_H[6] + eigvalue_c2*CnI_H[7] + (1/eigvalue_c3)*CnI_H[8] + CnI_H[9]; - const double L = std::abs(Utilities::fixed_power<2>(eigvalue_a1)*CnI_L[0] + eigvalue_a2*CnI_L[1] + (1/eigvalue_a3)*CnI_L[2] + Utilities::fixed_power<2>(eigvalue_b1)*CnI_L[3] + eigvalue_b2*CnI_L[4] + (1/eigvalue_b3)*CnI_L[5] + Utilities::fixed_power<2>(eigvalue_c1)*CnI_L[6] + eigvalue_c2*CnI_L[7] + (1/eigvalue_c3)*CnI_L[8] + CnI_L[9]); - const double M = std::abs(Utilities::fixed_power<2>(eigvalue_a1)*CnI_M[0] + eigvalue_a2*CnI_M[1] + (1/eigvalue_a3)*CnI_M[2] + Utilities::fixed_power<2>(eigvalue_b1)*CnI_M[3] + eigvalue_b2*CnI_M[4] + (1/eigvalue_b3)*CnI_M[5] + Utilities::fixed_power<2>(eigvalue_c1)*CnI_M[6] + eigvalue_c2*CnI_M[7] + (1/eigvalue_c3)*CnI_M[8] + CnI_M[9]); - const double N = std::abs(Utilities::fixed_power<2>(eigvalue_a1)*CnI_N[0] + eigvalue_a2*CnI_N[1] + (1/eigvalue_a3)*CnI_N[2] + Utilities::fixed_power<2>(eigvalue_b1)*CnI_N[3] + eigvalue_b2*CnI_N[4] + (1/eigvalue_b3)*CnI_N[5] + Utilities::fixed_power<2>(eigvalue_c1)*CnI_N[6] + eigvalue_c2*CnI_N[7] + (1/eigvalue_c3)*CnI_N[8] + CnI_N[9]); + // CPO2Hill v3 model: + const double F = Utilities::fixed_power<2>(eigvalue_a1)*CnI_F[0] + eigvalue_a1*CnI_F[1] + eigvalue_a2*CnI_F[2] + (1/eigvalue_a3)*CnI_F[3] + Utilities::fixed_power<2>(eigvalue_b1)*CnI_F[4] + eigvalue_b1*CnI_F[5] + eigvalue_b2*CnI_F[6] + (1/eigvalue_b3)*CnI_F[7] + Utilities::fixed_power<2>(eigvalue_c1)*CnI_F[8] + eigvalue_c1*CnI_F[9] + eigvalue_c2*CnI_F[10] + (1/eigvalue_c3)*CnI_F[11] + CnI_F[12]; + const double G = Utilities::fixed_power<2>(eigvalue_a1)*CnI_G[0] + eigvalue_a1*CnI_G[1] + eigvalue_a2*CnI_G[2] + (1/eigvalue_a3)*CnI_G[3] + Utilities::fixed_power<2>(eigvalue_b1)*CnI_G[4] + eigvalue_b1*CnI_G[5] + eigvalue_b2*CnI_G[6] + (1/eigvalue_b3)*CnI_G[7] + Utilities::fixed_power<2>(eigvalue_c1)*CnI_G[8] + eigvalue_c1*CnI_G[9] + eigvalue_c2*CnI_G[10] + (1/eigvalue_c3)*CnI_G[11] + CnI_G[12]; + const double H = Utilities::fixed_power<2>(eigvalue_a1)*CnI_H[0] + eigvalue_a1*CnI_H[1] + eigvalue_a2*CnI_H[2] + (1/eigvalue_a3)*CnI_H[3] + Utilities::fixed_power<2>(eigvalue_b1)*CnI_H[4] + eigvalue_b1*CnI_H[5] + eigvalue_b2*CnI_H[6] + (1/eigvalue_b3)*CnI_H[7] + Utilities::fixed_power<2>(eigvalue_c1)*CnI_H[8] + eigvalue_c1*CnI_H[9] + eigvalue_c2*CnI_H[10] + (1/eigvalue_c3)*CnI_H[11] + CnI_H[12]; + const double L = std::abs(Utilities::fixed_power<2>(eigvalue_a1)*CnI_L[0] + eigvalue_a1*CnI_L[1] + eigvalue_a2*CnI_L[2] + (1/eigvalue_a3)*CnI_L[3] + Utilities::fixed_power<2>(eigvalue_b1)*CnI_L[4] + eigvalue_b1*CnI_L[5] + eigvalue_b2*CnI_L[6] + (1/eigvalue_b3)*CnI_L[7] + Utilities::fixed_power<2>(eigvalue_c1)*CnI_L[8] + eigvalue_c1*CnI_L[9] + eigvalue_c2*CnI_L[10] + (1/eigvalue_c3)*CnI_L[11] + CnI_L[12]); + const double M = std::abs(Utilities::fixed_power<2>(eigvalue_a1)*CnI_M[0] + eigvalue_a1*CnI_M[1] + eigvalue_a2*CnI_M[2] + (1/eigvalue_a3)*CnI_M[3] + Utilities::fixed_power<2>(eigvalue_b1)*CnI_M[4] + eigvalue_b1*CnI_M[5] + eigvalue_b2*CnI_M[6] + (1/eigvalue_b3)*CnI_M[7] + Utilities::fixed_power<2>(eigvalue_c1)*CnI_M[8] + eigvalue_c1*CnI_M[9] + eigvalue_c2*CnI_M[10] + (1/eigvalue_c3)*CnI_M[11] + CnI_M[12]); + const double N = std::abs(Utilities::fixed_power<2>(eigvalue_a1)*CnI_N[0] + eigvalue_a1*CnI_N[1] + eigvalue_a2*CnI_N[2] + (1/eigvalue_a3)*CnI_N[3] + Utilities::fixed_power<2>(eigvalue_b1)*CnI_N[4] + eigvalue_b1*CnI_N[5] + eigvalue_b2*CnI_N[6] + (1/eigvalue_b3)*CnI_N[7] + Utilities::fixed_power<2>(eigvalue_c1)*CnI_N[8] + eigvalue_c1*CnI_N[9] + eigvalue_c2*CnI_N[10] + (1/eigvalue_c3)*CnI_N[11] + CnI_N[12]); Tensor<2,6> R_CPO_K; R_CPO_K[0][0] = Utilities::fixed_power<2>(R[0][0]); @@ -267,19 +267,19 @@ namespace aspect R_CPO_K[5][5] = R[0][0]*R[1][1]+R[0][1]*R[1][0]; SymmetricTensor<2,6> A; - A[0][0] = 2.0/3.0*(F+H); - A[0][1] = 2.0/3.0*(-F); - A[0][2] = 2.0/3.0*(-H); - A[1][1] = 2.0/3.0*(G+F); - A[1][2] = 2.0/3.0*(-G); - A[2][2] = 2.0/3.0*(H+G); - A[3][3] = 2.0/3.0*L; - A[4][4] = 2.0/3.0*M; - A[5][5] = 2.0/3.0*N; + A[0][0] = (G+H); + A[0][1] = (-H); + A[0][2] = (-G); + A[1][1] = (H+F); + A[1][2] = (-F); + A[2][2] = (F+G); + A[3][3] = (L); + A[4][4] = (M); + A[5][5] = (N); // A is the anisotropic tensor for the fluidity. We need its inverse, but it's not invertible due to singularity. // Thus we compute the Moore-Penrose pseudo inverse using SVD - LAPACKFullMatrix A_mat_lapack(6, 6), pinvA_mat_lapack(6,6); + LAPACKFullMatrix A_mat_lapack(6,6), pinvA_mat_lapack(6,6); for (unsigned int ai=0; ai<6; ++ai) { for (unsigned int aj=0; aj<6; ++aj) @@ -299,7 +299,7 @@ namespace aspect } // Calculate the fluidity tensor in the CPO frame - const Tensor<2,6> V = R_CPO_K * invA * transpose(R_CPO_K); + const Tensor<2,6> V = transpose(R_CPO_K) * invA * R_CPO_K; // Convert rank 2 viscosity tensor to rank 4 FullMatrix V_mat(6,6); @@ -333,17 +333,16 @@ namespace aspect // in which the stress is in MPa. SymmetricTensor<2,dim> stress = 2 * scalar_viscosity * V_r4 * deviatoric_strain_rate / 1e6; - const Tensor<2,dim> R_T = transpose(R); while (std::abs(residual) > threshold && n_iterations < max_iteration) { stress = (1./2.) * (stress + 2*scalar_viscosity * V_r4 * deviatoric_strain_rate / 1e6); - const Tensor<2,dim> S_CPO= R_T * stress * R; + const Tensor<2,dim> S_CPO= R * stress * transpose(R); - double Jhill = F*Utilities::fixed_power<2>(S_CPO[0][0]-S_CPO[1][1]) + G*Utilities::fixed_power<2>(S_CPO[1][1]-S_CPO[2][2]) + H*Utilities::fixed_power<2>(S_CPO[2][2]-S_CPO[0][0]) + 2*L*Utilities::fixed_power<2>(S_CPO[1][2]) + 2*M*Utilities::fixed_power<2>(S_CPO[0][2]) + 2*N*Utilities::fixed_power<2>(S_CPO[0][1]); + double Jhill = 2./3. * (F*Utilities::fixed_power<2>(S_CPO[1][1]-S_CPO[2][2]) + G*Utilities::fixed_power<2>(S_CPO[2][2]-S_CPO[0][0]) + H*Utilities::fixed_power<2>(S_CPO[0][0]-S_CPO[1][1]) + 2*L*Utilities::fixed_power<2>(S_CPO[1][2]) + 2*M*Utilities::fixed_power<2>(S_CPO[0][2]) + 2*N*Utilities::fixed_power<2>(S_CPO[0][1])); if (Jhill < 0) { - Jhill = std::abs(F)*Utilities::fixed_power<2>(S_CPO[0][0]-S_CPO[1][1]) + std::abs(G)*Utilities::fixed_power<2>(S_CPO[1][1]-S_CPO[2][2]) + std::abs(H)*Utilities::fixed_power<2>(S_CPO[2][2]-S_CPO[0][0]) + 2*L*Utilities::fixed_power<2>(S_CPO[1][2]) + 2*M*Utilities::fixed_power<2>(S_CPO[0][2]) + 2*N*Utilities::fixed_power<2>(S_CPO[0][1]); + Jhill = 2./3. * (std::abs(F)*Utilities::fixed_power<2>(S_CPO[1][1]-S_CPO[2][2]) + std::abs(G)*Utilities::fixed_power<2>(S_CPO[2][2]-S_CPO[0][0]) + std::abs(H)*Utilities::fixed_power<2>(S_CPO[0][0]-S_CPO[1][1]) + 2*L*Utilities::fixed_power<2>(S_CPO[1][2]) + 2*M*Utilities::fixed_power<2>(S_CPO[0][2]) + 2*N*Utilities::fixed_power<2>(S_CPO[0][1])); } AssertThrow(std::isfinite(Jhill), @@ -351,7 +350,7 @@ namespace aspect AssertThrow(Jhill >= 0, ExcMessage("Jhill should not be negative")); - const double scalar_viscosity_new = (1 / (Gamma * std::pow(Jhill,(n-1)/2))); + const double scalar_viscosity_new = (1 / ((2./3.) * Gamma * std::pow(Jhill,(n-1)/2))); residual = std::abs(scalar_viscosity_new - scalar_viscosity); scalar_viscosity = scalar_viscosity_new; threshold = 0.001*scalar_viscosity; @@ -373,15 +372,15 @@ namespace aspect { // Assign an isotropic viscosity tensor SymmetricTensor<2,6> V; - V[0][0] = 2.0/3.0; - V[0][1] = -1.0/3.0; - V[0][2] = -1.0/3.0; - V[1][1] = 2.0/3.0; - V[1][2] = -1.0/3.0; + V[0][0] = 4.0/9.0; + V[0][1] = -2.0/9.0; + V[0][2] = -2.0/9.0; + V[1][1] = 4.0/9.0; + V[1][2] = -2.0/9.0; V[2][2] = 2.0/3.0; - V[3][3] = 1; - V[4][4] = 1; - V[5][5] = 1; + V[3][3] = 2.0/3.0; + V[4][4] = 2.0/3.0; + V[5][5] = 2.0/3.0; // Convert rank 2 viscosity tensor to rank 4 FullMatrix V_mat(6,6); @@ -428,9 +427,9 @@ namespace aspect { EquationOfState::LinearizedIncompressible::declare_parameters (prm); - prm.declare_entry ("Coefficients and intercept for F", "1.0390459583037057, -0.767458622, 0.003066208, 0.19651133418307049, 0.413093763, 0.015463162, -0.935925291, -2.392877563, 0.051834768, 1.0799807050187482", + prm.declare_entry ("Coefficients and intercept for F", "0.5920219168461529, -0.831936049, -0.000937583, -0.00029648, 0.380413345, -0.533048795, 0.46835862365145753, -0.000965503, -1.249340274, 1.0748554477472438, -0.167662132, 0.003358407, 0.5215020195972386", Patterns::List(Patterns::Double()), - "9 Coefficients and 1 intercept to compute the Hill Parameter F " + "12 Coefficients and 1 intercept to compute the Hill Parameter F " "according to the linear regression relation provided in the cookbook documentation. " "The first 3 coefficients are multiplied respectively by: " "the square of the largest eigenvalue, the second-largest eigenvalue, " @@ -438,21 +437,21 @@ namespace aspect "The next 3 coefficients are used in the same way for the b-axis, " "and the final 3 for the c-axis. Together with the intercept, " "these values form the full regression expression for F."); - prm.declare_entry ("Coefficients and intercept for G", "-2.836270315, -1.632453092, 0.000687606, 0.2671850239576621, -0.993392913, 0.002699241, 1.9689530759060374, 2.314442451425019, -0.018655905, 0.6887411607403755", + prm.declare_entry ("Coefficients and intercept for G", "-1.6951323, 1.3364976547800977, -0.18410694, 4.918308228973878e-05, 0.7501414478371807, 0.691412915, 0.37696216069289673, -0.001537058, -0.66969478, -0.551507796, -0.428462988, 0.003403174, 0.2602865312446454", Patterns::List(Patterns::Double()), - "9 Coefficients and 1 intercept to compute the Hill Parameter G in the same way as above."); - prm.declare_entry ("Coefficients and intercept for H", "1.6687493021559732, 0.5797579293682223, 0.003241593, 0.701661336, 0.2513824481429968, 0.000229291, -2.003227619, -2.57032429, 0.071454541, 0.7490268673620638", + "12 Coefficients and 1 intercept to compute the Hill Parameter G in the same way as above."); + prm.declare_entry ("Coefficients and intercept for H", "-1.139612048, 1.353113344145978, 0.7510486623018213, -0.001656848, -0.255721452, -1.006433455, -0.11595106, 0.003177149, 0.6837240306536184, -0.031162568, -0.080356281, 0.005621241, 0.26788414888354445", Patterns::List(Patterns::Double()), - "9 Coefficients and 1 intercept to compute the Hill Parameter H in the same way as above."); - prm.declare_entry ("Coefficients and intercept for L", "-0.325145943, 0.7284642859944138, 0.000404879, -0.665446098, 0.5152847961409479, 0.002722782, -1.026786493, -1.262574542, 0.009168498, 1.595422603", + "12 Coefficients and 1 intercept to compute the Hill Parameter H in the same way as above."); + prm.declare_entry ("Coefficients and intercept for L", "-3.510950516, 2.6864808033885543, 0.035838123, -0.000504338, 3.9483598066088383, -3.816102334, -0.778569714, 0.003688104, 4.122460734346824, -2.482527095, 1.3200590504614733, -0.002399896, 2.0027068994912076", Patterns::List(Patterns::Double()), - "9 Coefficients and 1 intercept to compute the Hill Parameter L in the same way as above."); - prm.declare_entry ("Coefficients and intercept for M", "1.6427437063774875, 0.8777500120437522, 0.004651732, 2.489417876177839, 0.8162729707609052, -0.010736521, -2.49420455, -0.511446494, -0.009362491, 0.893677343", + "12 Coefficients and 1 intercept to compute the Hill Parameter L in the same way as above."); + prm.declare_entry ("Coefficients and intercept for M", "4.536980494567378, -3.227568914, 0.27609495132676254, 0.007436169, -7.446913908, 5.763821882498737, -1.4026181, 0.032132134, 2.9678024468288697, -3.434721081, -2.265560577, 0.1215179917888699, 2.4409386788998457", Patterns::List(Patterns::Double()), - "9 Coefficients and 1 intercept to compute the Hill Parameter M in the same way as above."); - prm.declare_entry ("Coefficients and intercept for N", "0.8122098589701904, 0.15663795996228266, 0.001500252, -1.648578168, 0.19362392490527092, -0.009650519, 1.6796559729985163, -0.103640482, 0.01971017, 1.2132200780065174", + "12 Coefficients and 1 intercept to compute the Hill Parameter M in the same way as above."); + prm.declare_entry ("Coefficients and intercept for N", "7.872922986831338, -7.933513948, -2.588175191, 0.029843804645040883, 7.605864624755694, -5.469451775, -0.347688637, 0.06395131, -1.78763278, 2.2550636824173584, 3.023166891521831, -0.102862765, 3.6958741003498234", Patterns::List(Patterns::Double()), - "9 Coefficients and 1 intercept to compute the Hill Parameter N in the same way as above."); + "12 Coefficients and 1 intercept to compute the Hill Parameter N in the same way as above."); prm.declare_entry ("Reference viscosity", "1e9", Patterns::Double(), diff --git a/cookbooks/CPO_induced_anisotropic_viscosity/shearbox_cpo_av.prm b/cookbooks/CPO_induced_anisotropic_viscosity/shearbox_cpo_av.prm index e2fd3a6561e..e68ec940fc3 100644 --- a/cookbooks/CPO_induced_anisotropic_viscosity/shearbox_cpo_av.prm +++ b/cookbooks/CPO_induced_anisotropic_viscosity/shearbox_cpo_av.prm @@ -10,7 +10,7 @@ set Dimension = 3 set Pressure normalization = surface set Surface pressure = 0 set Nonlinear solver scheme = single Advection, single Stokes -set End time = 0.1 +set End time = 0.5 set Use years instead of seconds = false set Output directory = output_Shearbox_CPO_AV @@ -119,7 +119,7 @@ subsection Postprocess subsection Visualization set Time between graphical output = 0.1 - set List of output variables = material properties, strain rate, named additional outputs, shear stress, stress + set List of output variables = material properties, strain rate, named additional outputs, shear stress, stress, Anisotropic stress subsection Material properties set List of material properties = density, viscosity diff --git a/doc/sphinx/references.bib b/doc/sphinx/references.bib index ef459de1f04..7e8ebed47e9 100644 --- a/doc/sphinx/references.bib +++ b/doc/sphinx/references.bib @@ -12671,4 +12671,19 @@ @article{bingham:1974 number={6}, pages={1201--1225}, year={1974} -} \ No newline at end of file +} + +@software{kiraly:etal:2026, + author = {Agi Kiraly and + Yijun Wang and + Clinton P. Conrad and + Lars Hansen and + Ben Mather}, + title = {Dataset for the paper titled 'CPO2Hill: A new model linking olivine texture parameters to anisotropic viscous behavior' by Király et al.}, + month = jan, + year = 2026, + publisher = {Zenodo}, + version = {v2}, + doi = {10.5281/zenodo.18611490}, + url = {https://doi.org/10.5281/zenodo.3924604} +} \ No newline at end of file diff --git a/tests/cpo_induced_anisotropic_viscosity.prm b/tests/cpo_induced_anisotropic_viscosity.prm index f49e41d4dc7..d7aa4140132 100644 --- a/tests/cpo_induced_anisotropic_viscosity.prm +++ b/tests/cpo_induced_anisotropic_viscosity.prm @@ -1,7 +1,7 @@ include $ASPECT_SOURCE_DIR/cookbooks/CPO_induced_anisotropic_viscosity/shearbox_cpo_av.prm set Additional shared libraries = ./libCPO_induced_anisotropic_viscosity.debug.so -set End time = 0 +set End time = 0.5 subsection Particles subsection Generator diff --git a/tests/cpo_induced_anisotropic_viscosity/screen-output b/tests/cpo_induced_anisotropic_viscosity/screen-output index 526bda4233f..493da49f8ec 100644 --- a/tests/cpo_induced_anisotropic_viscosity/screen-output +++ b/tests/cpo_induced_anisotropic_viscosity/screen-output @@ -1,5 +1,5 @@ -Loading shared library <./libcpo_induced_anisotropic_viscosity.debug.so> +Loading shared library <./libcpo_induced_anisotropic_viscosity.release.so> Number of active cells: 1 (on 1 levels) Number of degrees of freedom: 494 (81+8+27+27+27+27+27+27+27+27+27+27+27+27+27+27+27) @@ -9,7 +9,7 @@ Number of degrees of freedom: 494 (81+8+27+27+27+27+27+27+27+27+27+27+27+27+27+2 Advecting particles... done. Copying properties into prescribed compositional field scalar_viscosity... done. Rebuilding Stokes preconditioner... - Solving Stokes system (AMG)... 3+0 iterations. + Solving Stokes system (AMG)... 4+0 iterations. Postprocessing: RMS, max velocity: 0.289 m/s, 0.5 m/s @@ -18,6 +18,76 @@ Number of degrees of freedom: 494 (81+8+27+27+27+27+27+27+27+27+27+27+27+27+27+2 Writing particle output: output-cpo_induced_anisotropic_viscosity/particles/particles-00000 Writing particle cpo output: output-cpo_induced_anisotropic_viscosity/particles_cpo/CPO-00000 +*** Timestep 1: t=0.1 seconds, dt=0.1 seconds + Solving temperature system... 0 iterations. + Advecting particles... done. + Copying properties into prescribed compositional field scalar_viscosity... done. + Rebuilding Stokes preconditioner... + Solving Stokes system (AMG)... 9+0 iterations. + + Postprocessing: + RMS, max velocity: 0.289 m/s, 0.5 m/s + Compositions min/max/mass: 1.864e+07/1.864e+07/1.864e+07 // 4.654/4.654/4.654 // 0.6551/0.6551/0.6551 // 0.2233/0.2233/0.2233 // 0.1217/0.1217/0.1217 // 2.055/2.055/2.055 // 0.5386/0.5386/0.5386 // 0.4457/0.4457/0.4457 // 0.01567/0.01567/0.01567 // 0.08621/0.08621/0.08621 // 0.4653/0.4653/0.4653 // 0.3568/0.3568/0.3568 // 0.1779/0.1779/0.1779 // 0/0/0 + Writing graphical output: output-cpo_induced_anisotropic_viscosity/solution/solution-00001 + Writing particle output: output-cpo_induced_anisotropic_viscosity/particles/particles-00001 + Writing particle cpo output: output-cpo_induced_anisotropic_viscosity/particles_cpo/CPO-00001 + +*** Timestep 2: t=0.2 seconds, dt=0.1 seconds + Solving temperature system... 0 iterations. + Advecting particles... done. + Copying properties into prescribed compositional field scalar_viscosity... done. + Rebuilding Stokes preconditioner... + Solving Stokes system (AMG)... 13+0 iterations. + + Postprocessing: + RMS, max velocity: 0.289 m/s, 0.5 m/s + Compositions min/max/mass: 1.635e+07/1.635e+07/1.635e+07 // 4.682/4.682/4.682 // 0.687/0.687/0.687 // 0.2543/0.2543/0.2543 // 0.05866/0.05866/0.05866 // 1.956/1.956/1.956 // 0.6427/0.6427/0.6427 // 0.318/0.318/0.318 // 0.03931/0.03931/0.03931 // 0.2642/0.2642/0.2642 // 0.5393/0.5393/0.5393 // 0.2979/0.2979/0.2979 // 0.1628/0.1628/0.1628 // 0/0/0 + Writing graphical output: output-cpo_induced_anisotropic_viscosity/solution/solution-00002 + Writing particle output: output-cpo_induced_anisotropic_viscosity/particles/particles-00002 + Writing particle cpo output: output-cpo_induced_anisotropic_viscosity/particles_cpo/CPO-00002 + +*** Timestep 3: t=0.3 seconds, dt=0.1 seconds + Solving temperature system... 0 iterations. + Advecting particles... done. + Copying properties into prescribed compositional field scalar_viscosity... done. + Rebuilding Stokes preconditioner... + Solving Stokes system (AMG)... 7+0 iterations. + + Postprocessing: + RMS, max velocity: 0.289 m/s, 0.5 m/s + Compositions min/max/mass: 1.636e+07/1.637e+07/1.637e+07 // 4.901/4.901/4.901 // 0.5217/0.5217/0.5217 // 0.3452/0.3452/0.3452 // 0.1331/0.1331/0.1331 // 1.327/1.327/1.327 // 0.507/0.507/0.507 // 0.4874/0.4874/0.4874 // 0.005555/0.005555/0.005555 // 6.167/6.167/6.167 // 0.8138/0.8138/0.8138 // 0.165/0.165/0.165 // 0.02117/0.02117/0.02117 // 0/0/0 + Writing graphical output: output-cpo_induced_anisotropic_viscosity/solution/solution-00003 + Writing particle output: output-cpo_induced_anisotropic_viscosity/particles/particles-00003 + Writing particle cpo output: output-cpo_induced_anisotropic_viscosity/particles_cpo/CPO-00003 + +*** Timestep 4: t=0.4 seconds, dt=0.1 seconds + Solving temperature system... 0 iterations. + Advecting particles... done. + Copying properties into prescribed compositional field scalar_viscosity... done. + Rebuilding Stokes preconditioner... + Solving Stokes system (AMG)... 6+0 iterations. + + Postprocessing: + RMS, max velocity: 0.289 m/s, 0.5 m/s + Compositions min/max/mass: 1.814e+07/1.814e+07/1.814e+07 // 5.401/5.401/5.401 // 0.5239/0.5239/0.5239 // 0.3726/0.3726/0.3726 // 0.1035/0.1035/0.1035 // 2.101/2.101/2.101 // 0.5855/0.5855/0.5855 // 0.3469/0.3469/0.3469 // 0.06767/0.06767/0.06767 // 4.669/4.669/4.669 // 0.674/0.674/0.674 // 0.1845/0.1845/0.1845 // 0.1415/0.1415/0.1415 // 0/0/0 + Writing graphical output: output-cpo_induced_anisotropic_viscosity/solution/solution-00004 + Writing particle output: output-cpo_induced_anisotropic_viscosity/particles/particles-00004 + Writing particle cpo output: output-cpo_induced_anisotropic_viscosity/particles_cpo/CPO-00004 + +*** Timestep 5: t=0.5 seconds, dt=0.1 seconds + Solving temperature system... 0 iterations. + Advecting particles... done. + Copying properties into prescribed compositional field scalar_viscosity... done. + Rebuilding Stokes preconditioner... + Solving Stokes system (AMG)... 6+0 iterations. + + Postprocessing: + RMS, max velocity: 0.289 m/s, 0.5 m/s + Compositions min/max/mass: 1.758e+07/1.758e+07/1.758e+07 // 4.104/4.104/4.104 // 0.7394/0.7394/0.7394 // 0.1602/0.1602/0.1602 // 0.1004/0.1004/0.1004 // 1.107/1.107/1.107 // 0.5986/0.5986/0.5986 // 0.3103/0.3103/0.3103 // 0.09108/0.09108/0.09108 // 6.233/6.233/6.233 // 0.563/0.563/0.563 // 0.2924/0.2924/0.2924 // 0.1446/0.1446/0.1446 // 0/0/0 + Writing graphical output: output-cpo_induced_anisotropic_viscosity/solution/solution-00005 + Writing particle output: output-cpo_induced_anisotropic_viscosity/particles/particles-00005 + Writing particle cpo output: output-cpo_induced_anisotropic_viscosity/particles_cpo/CPO-00005 + Termination requested by criterion: end time diff --git a/tests/cpo_induced_anisotropic_viscosity/statistics b/tests/cpo_induced_anisotropic_viscosity/statistics index 2082c7e4b56..8fdacf89586 100644 --- a/tests/cpo_induced_anisotropic_viscosity/statistics +++ b/tests/cpo_induced_anisotropic_viscosity/statistics @@ -59,4 +59,9 @@ # 59: Number of advected particles # 60: Particle file name # 61: Particle CPO file name -0 0.000000000000e+00 0.000000000000e+00 1 89 27 378 0 0 0 3 3 3 2.88675135e-01 5.00000000e-01 1.00000000e+09 1.00000000e+09 1.00000000e+09 5.27630079e+00 5.27630079e+00 5.27630079e+00 5.20309402e-01 5.20309402e-01 5.20309402e-01 3.47923711e-01 3.47923711e-01 3.47923711e-01 1.31766887e-01 1.31766887e-01 1.31766887e-01 1.73622715e+00 1.73622715e+00 1.73622715e+00 5.74325744e-01 5.74325744e-01 5.74325744e-01 3.80780388e-01 3.80780388e-01 3.80780388e-01 4.48938682e-02 4.48938682e-02 4.48938682e-02 5.43504380e+00 5.43504380e+00 5.43504380e+00 6.07858069e-01 6.07858069e-01 6.07858069e-01 3.10569263e-01 3.10569263e-01 3.10569263e-01 8.15726682e-02 8.15726682e-02 8.15726682e-02 0.00000000e+00 0.00000000e+00 0.00000000e+00 output-cpo_induced_anisotropic_viscosity/solution/solution-00000 1 output-cpo_induced_anisotropic_viscosity/particles/particles-00000 output-cpo_induced_anisotropic_viscosity/particles_cpo/CPO-00000 +0 0.000000000000e+00 0.000000000000e+00 1 89 27 378 0 0 0 4 4 4 2.88675135e-01 5.00000000e-01 1.00000000e+09 1.00000000e+09 1.00000000e+09 5.27630079e+00 5.27630079e+00 5.27630079e+00 5.20309402e-01 5.20309402e-01 5.20309402e-01 3.47923711e-01 3.47923711e-01 3.47923711e-01 1.31766887e-01 1.31766887e-01 1.31766887e-01 1.73622715e+00 1.73622715e+00 1.73622715e+00 5.74325744e-01 5.74325744e-01 5.74325744e-01 3.80780388e-01 3.80780388e-01 3.80780388e-01 4.48938682e-02 4.48938682e-02 4.48938682e-02 5.43504380e+00 5.43504380e+00 5.43504380e+00 6.07858069e-01 6.07858069e-01 6.07858069e-01 3.10569263e-01 3.10569263e-01 3.10569263e-01 8.15726682e-02 8.15726682e-02 8.15726682e-02 0.00000000e+00 0.00000000e+00 0.00000000e+00 output-cpo_induced_anisotropic_viscosity/solution/solution-00000 1 output-cpo_induced_anisotropic_viscosity/particles/particles-00000 output-cpo_induced_anisotropic_viscosity/particles_cpo/CPO-00000 +1 1.000000000000e-01 1.000000000000e-01 1 89 27 378 0 0 0 9 9 9 2.88675135e-01 5.00000000e-01 1.86352107e+07 1.86352107e+07 1.86352107e+07 4.65379362e+00 4.65379362e+00 4.65379362e+00 6.55064477e-01 6.55064477e-01 6.55064477e-01 2.23254860e-01 2.23254860e-01 2.23254860e-01 1.21680662e-01 1.21680662e-01 1.21680662e-01 2.05523620e+00 2.05523620e+00 2.05523620e+00 5.38639659e-01 5.38639659e-01 5.38639659e-01 4.45692291e-01 4.45692291e-01 4.45692291e-01 1.56680496e-02 1.56680496e-02 1.56680496e-02 8.62077695e-02 8.62077695e-02 8.62077695e-02 4.65288335e-01 4.65288335e-01 4.65288335e-01 3.56763907e-01 3.56763907e-01 3.56763907e-01 1.77947758e-01 1.77947758e-01 1.77947758e-01 0.00000000e+00 0.00000000e+00 0.00000000e+00 output-cpo_induced_anisotropic_viscosity/solution/solution-00001 1 output-cpo_induced_anisotropic_viscosity/particles/particles-00001 output-cpo_induced_anisotropic_viscosity/particles_cpo/CPO-00001 +2 2.000000000000e-01 1.000000000000e-01 1 89 27 378 0 0 0 13 13 13 2.88675135e-01 5.00000000e-01 1.63538553e+07 1.63538557e+07 1.63538555e+07 4.68191169e+00 4.68191169e+00 4.68191169e+00 6.86995446e-01 6.86995446e-01 6.86995446e-01 2.54340432e-01 2.54340432e-01 2.54340432e-01 5.86641214e-02 5.86641214e-02 5.86641214e-02 1.95593521e+00 1.95593521e+00 1.95593521e+00 6.42677563e-01 6.42677563e-01 6.42677563e-01 3.18014417e-01 3.18014417e-01 3.18014417e-01 3.93080207e-02 3.93080207e-02 3.93080207e-02 2.64223298e-01 2.64223298e-01 2.64223298e-01 5.39265192e-01 5.39265192e-01 5.39265192e-01 2.97927889e-01 2.97927889e-01 2.97927889e-01 1.62806919e-01 1.62806919e-01 1.62806919e-01 0.00000000e+00 0.00000000e+00 0.00000000e+00 output-cpo_induced_anisotropic_viscosity/solution/solution-00002 1 output-cpo_induced_anisotropic_viscosity/particles/particles-00002 output-cpo_induced_anisotropic_viscosity/particles_cpo/CPO-00002 +3 3.000000000000e-01 1.000000000000e-01 1 89 27 378 0 0 0 7 7 7 2.88675135e-01 5.00000000e-01 1.63648179e+07 1.63660690e+07 1.63654434e+07 4.90124042e+00 4.90124042e+00 4.90124042e+00 5.21650669e-01 5.21650669e-01 5.21650669e-01 3.45223725e-01 3.45223725e-01 3.45223725e-01 1.33125606e-01 1.33125606e-01 1.33125606e-01 1.32703227e+00 1.32703227e+00 1.32703227e+00 5.07011584e-01 5.07011584e-01 5.07011584e-01 4.87433154e-01 4.87433154e-01 4.87433154e-01 5.55526238e-03 5.55526238e-03 5.55526238e-03 6.16742540e+00 6.16742540e+00 6.16742540e+00 8.13791350e-01 8.13791350e-01 8.13791350e-01 1.65040519e-01 1.65040519e-01 1.65040519e-01 2.11681311e-02 2.11681311e-02 2.11681311e-02 0.00000000e+00 0.00000000e+00 0.00000000e+00 output-cpo_induced_anisotropic_viscosity/solution/solution-00003 1 output-cpo_induced_anisotropic_viscosity/particles/particles-00003 output-cpo_induced_anisotropic_viscosity/particles_cpo/CPO-00003 +4 4.000000000000e-01 1.000000000000e-01 1 89 27 378 0 0 0 6 6 6 2.88675135e-01 5.00000000e-01 1.81430453e+07 1.81430495e+07 1.81430474e+07 5.40122005e+00 5.40122005e+00 5.40122005e+00 5.23943612e-01 5.23943612e-01 5.23943612e-01 3.72586477e-01 3.72586477e-01 3.72586477e-01 1.03469912e-01 1.03469912e-01 1.03469912e-01 2.10070182e+00 2.10070182e+00 2.10070182e+00 5.85454430e-01 5.85454430e-01 5.85454430e-01 3.46873102e-01 3.46873102e-01 3.46873102e-01 6.76724685e-02 6.76724685e-02 6.76724685e-02 4.66879242e+00 4.66879242e+00 4.66879242e+00 6.73983116e-01 6.73983116e-01 6.73983116e-01 1.84471150e-01 1.84471150e-01 1.84471150e-01 1.41545734e-01 1.41545734e-01 1.41545734e-01 0.00000000e+00 0.00000000e+00 0.00000000e+00 output-cpo_induced_anisotropic_viscosity/solution/solution-00004 1 output-cpo_induced_anisotropic_viscosity/particles/particles-00004 output-cpo_induced_anisotropic_viscosity/particles_cpo/CPO-00004 +5 5.000000000000e-01 1.000000000000e-01 1 89 27 378 0 0 0 6 6 6 2.88675135e-01 5.00000000e-01 1.75833622e+07 1.75833623e+07 1.75833622e+07 4.10371973e+00 4.10371973e+00 4.10371973e+00 7.39363795e-01 7.39363795e-01 7.39363795e-01 1.60220058e-01 1.60220058e-01 1.60220058e-01 1.00416147e-01 1.00416147e-01 1.00416147e-01 1.10714450e+00 1.10714450e+00 1.10714450e+00 5.98598415e-01 5.98598415e-01 5.98598415e-01 3.10324521e-01 3.10324521e-01 3.10324521e-01 9.10770642e-02 9.10770642e-02 9.10770642e-02 6.23337747e+00 6.23337747e+00 6.23337747e+00 5.63042860e-01 5.63042860e-01 5.63042860e-01 2.92401912e-01 2.92401912e-01 2.92401912e-01 1.44555227e-01 1.44555227e-01 1.44555227e-01 0.00000000e+00 0.00000000e+00 0.00000000e+00 output-cpo_induced_anisotropic_viscosity/solution/solution-00005 1 output-cpo_induced_anisotropic_viscosity/particles/particles-00005 output-cpo_induced_anisotropic_viscosity/particles_cpo/CPO-00005