Skip to content

Commit

Permalink
Update test_christoffel_stepbystep.ipynb
Browse files Browse the repository at this point in the history
  • Loading branch information
@marcoalopez
marcoalopez committed Jul 3, 2024
1 parent a47d4a7 commit 0eb0b4f
Showing 1 changed file with 64 additions and 53 deletions.
117 changes: 64 additions & 53 deletions src/test_christoffel_stepbystep.ipynb
View file
Original file line number Diff line number Diff line change
Expand Up @@ -29,11 +29,13 @@
"source": [
"The goal of PyRockWave's Christoffel module is to estimate various seismic properties as a function of direction using the Christoffel equation (Christoffel, 1877). These properties include the $V_p$ (compressional), $V_{s1}$, $V_{s2}$ (shear) wave velocities (phase and group) in km/s, the polarisation directions of the sound waves, the shear wave splitting, the coordinates or the ray surface, the power flow angle and the enhancement factor. We will define all these parameters throughout this notebook.\n",
"\n",
"> _Other material properties unrelated to seismic waves are Young modulus, shear modulus, Poisson's ratio, linear compresibility and Debye temperature._\n",
"\n",
"## Input data\n",
"\n",
"The necessary input data to estimate the seismic properties in any direction of a given material using the Christoffel equation are the density and the stiffness tensor of the material.\n",
"The necessary input data to estimate the seismic and material properties in any direction of a given material using the Christoffel equation are the density and the stiffness tensor of the material.\n",
"\n",
"The stiffness tensor $C$ is a fundamental property of a material that generalises Hooke's law in three dimensions, relating strains to stresses in the elastic regime. This is usually given (or abbreviated) in a 6x6 matrix ($C_{ij}$) using Voigt's notation, and values are in $GPa$. Density is usually given in...TODO "
"The stiffness tensor $C$ is a fundamental property of a material that generalises Hooke's law in three dimensions, relating strains to stresses in the elastic regime. This is usually given (or abbreviated) in a 6x6 matrix ($C_{ij}$) using Voigt's notation, and values are in $GPa$. Density is usually given in $g/cm^3$."
]
},
{
Expand Down Expand Up @@ -61,7 +63,7 @@
"\n",
"For most materials, the seismic properties are anisotropic, i.e. they depend on the direction. In such cases, the Christoffel equation must be solved for any propagation direction of interest, defined as wave vectors or $q$...\n",
"\n",
"Wave vectors can be provides using spherical coordinates...TODO"
"Wave vectors can be defined using spherical coordinates (polar and azimuthal angles) or in Cartesian coordinates ($x, y, z$) normalized to unit sphere (i.e., vector of size 1). For calculations Cartesian coordinates will be used."
]
},
{
Expand All @@ -70,7 +72,7 @@
"metadata": {},
"outputs": [],
"source": [
"# create wavevectors in radians\n",
"# create equispaced wavevectors in spherical coordinates (in radians)\n",
"azimuths, polar = c.equispaced_S2_grid(n=150)"
]
},
Expand Down Expand Up @@ -118,7 +120,7 @@
}
],
"source": [
"# Create the wavevectors\n",
"# Store wavevectors in a single Numpy array for calculations\n",
"wavevectors = np.column_stack((x, y, z))\n",
"\n",
"wavevectors.shape"
Expand Down Expand Up @@ -249,19 +251,7 @@
"source": [
"## Step 4: Estimate the eigenvalues $\\lambda$ and eigenvectors\n",
"\n",
"The eigenvalues $\\lambda_1$, $\\lambda_2$ and $\\lambda_3$ of the Christoffel tensor $M_{ij}$ are always positive and related to the wave velocities $V_p$, $V_{s1}$, $V_{s2}$ propagating in the direction $\\vec{n}$ by the formulae\n",
"\n",
"$$\n",
"V_p = \\sqrt{ \\frac{\\lambda_1}{\\rho} }, \\quad V_{s1} = \\sqrt{ \\frac{\\lambda_2}{\\rho} }, \\quad V_{s2} = \\sqrt{ \\frac{\\lambda_3}{\\rho} }\n",
"$$\n",
"\n",
"where $\\rho$ denotes the material density. Note that if the Christoffel tensor is already normalized to the density this is just\n",
"\n",
"$$\n",
"v_p = \\sqrt{\\lambda}, \\quad \\lambda = v^2_p\n",
"$$\n",
"\n",
"where $v_p$ denotes phase velocities."
"TODO"
]
},
{
Expand Down Expand Up @@ -290,7 +280,21 @@
"cell_type": "markdown",
"metadata": {},
"source": [
"## Step 5: Calculate phase velocities ($v_p$)"
"## Step 5: Calculate phase velocities ($v_p$)\n",
"\n",
"The eigenvalues $\\lambda_1$, $\\lambda_2$ and $\\lambda_3$ of the Christoffel tensor $M_{ij}$ are always real, positive and related to the wave velocities $V_p$, $V_{s1}$, $V_{s2}$ propagating in the direction $\\vec{n}$ by the formulae\n",
"\n",
"$$\n",
"V_p = \\sqrt{ \\frac{\\lambda_1}{\\rho} }, \\quad V_{s1} = \\sqrt{ \\frac{\\lambda_2}{\\rho} }, \\quad V_{s2} = \\sqrt{ \\frac{\\lambda_3}{\\rho} }\n",
"$$\n",
"\n",
"where $\\rho$ denotes the material density. Note that if the Christoffel tensor is already normalized to the density this is just\n",
"\n",
"$$\n",
"v_p = \\sqrt{\\lambda_n}, \\quad \\lambda_n = v^2_p\n",
"$$\n",
"\n",
"where $v_p$ denotes phase velocities. TODO -> define phase velocities"
]
},
{
Expand Down Expand Up @@ -488,6 +492,13 @@
"group_velocity_magnitudes[:10]"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"> This wave velocities are not normalized by the density of the material!"
]
},
{
"cell_type": "markdown",
"metadata": {},
Expand Down Expand Up @@ -634,8 +645,8 @@
"outputs": [],
"source": [
"data = pd.DataFrame({\n",
" 'polar_ang': np.rad2deg(polar),\n",
" 'azimuthal_ang': np.rad2deg(azimuths),\n",
" 'polar_ang': np.around(np.rad2deg(polar), 1),\n",
" 'azimuthal_ang': np.around(np.rad2deg(azimuths), 1),\n",
" 'phase_Vp': phase_velocities[:, 2],\n",
" 'phase_Vs1': phase_velocities[:, 1],\n",
" 'phase_Vs2': phase_velocities[:, 0],\n",
Expand Down Expand Up @@ -684,8 +695,8 @@
" <tbody>\n",
" <tr>\n",
" <th>0</th>\n",
" <td>0.000000</td>\n",
" <td>0.000000</td>\n",
" <td>0.0</td>\n",
" <td>0.0</td>\n",
" <td>4.731335</td>\n",
" <td>3.345559</td>\n",
" <td>3.345559</td>\n",
Expand All @@ -695,8 +706,8 @@
" </tr>\n",
" <tr>\n",
" <th>1</th>\n",
" <td>16.930725</td>\n",
" <td>0.000000</td>\n",
" <td>16.9</td>\n",
" <td>0.0</td>\n",
" <td>4.932846</td>\n",
" <td>3.345559</td>\n",
" <td>3.040613</td>\n",
Expand All @@ -706,8 +717,8 @@
" </tr>\n",
" <tr>\n",
" <th>2</th>\n",
" <td>19.327569</td>\n",
" <td>222.492236</td>\n",
" <td>19.3</td>\n",
" <td>222.5</td>\n",
" <td>4.988421</td>\n",
" <td>3.262362</td>\n",
" <td>3.040346</td>\n",
Expand All @@ -717,8 +728,8 @@
" </tr>\n",
" <tr>\n",
" <th>3</th>\n",
" <td>21.467383</td>\n",
" <td>84.984472</td>\n",
" <td>21.5</td>\n",
" <td>85.0</td>\n",
" <td>5.013496</td>\n",
" <td>3.343000</td>\n",
" <td>2.908654</td>\n",
Expand All @@ -728,8 +739,8 @@
" </tr>\n",
" <tr>\n",
" <th>4</th>\n",
" <td>23.420812</td>\n",
" <td>307.476708</td>\n",
" <td>23.4</td>\n",
" <td>307.5</td>\n",
" <td>5.067123</td>\n",
" <td>3.234312</td>\n",
" <td>2.938460</td>\n",
Expand All @@ -750,8 +761,8 @@
" </tr>\n",
" <tr>\n",
" <th>145</th>\n",
" <td>156.579188</td>\n",
" <td>358.881977</td>\n",
" <td>156.6</td>\n",
" <td>358.9</td>\n",
" <td>5.045805</td>\n",
" <td>3.345408</td>\n",
" <td>2.849412</td>\n",
Expand All @@ -761,8 +772,8 @@
" </tr>\n",
" <tr>\n",
" <th>146</th>\n",
" <td>158.532617</td>\n",
" <td>221.374213</td>\n",
" <td>158.5</td>\n",
" <td>221.4</td>\n",
" <td>5.030398</td>\n",
" <td>3.244671</td>\n",
" <td>2.989694</td>\n",
Expand All @@ -772,8 +783,8 @@
" </tr>\n",
" <tr>\n",
" <th>147</th>\n",
" <td>160.672431</td>\n",
" <td>83.866449</td>\n",
" <td>160.7</td>\n",
" <td>83.9</td>\n",
" <td>4.976185</td>\n",
" <td>3.342438</td>\n",
" <td>2.972667</td>\n",
Expand All @@ -783,8 +794,8 @@
" </tr>\n",
" <tr>\n",
" <th>148</th>\n",
" <td>163.069275</td>\n",
" <td>306.358685</td>\n",
" <td>163.1</td>\n",
" <td>306.4</td>\n",
" <td>4.940651</td>\n",
" <td>3.287877</td>\n",
" <td>3.090450</td>\n",
Expand All @@ -794,8 +805,8 @@
" </tr>\n",
" <tr>\n",
" <th>149</th>\n",
" <td>180.000000</td>\n",
" <td>0.000000</td>\n",
" <td>180.0</td>\n",
" <td>0.0</td>\n",
" <td>4.731335</td>\n",
" <td>3.345559</td>\n",
" <td>3.345559</td>\n",
Expand All @@ -809,18 +820,18 @@
"</div>"
],
"text/plain": [
" polar_ang azimuthal_ang phase_Vp phase_Vs1 phase_Vs2 group_Vp \\\n",
"0 0.000000 0.000000 4.731335 3.345559 3.345559 25.109194 \n",
"1 16.930725 0.000000 4.932846 3.345559 3.040613 18.736040 \n",
"2 19.327569 222.492236 4.988421 3.262362 3.040346 25.801391 \n",
"3 21.467383 84.984472 5.013496 3.343000 2.908654 21.006679 \n",
"4 23.420812 307.476708 5.067123 3.234312 2.938460 16.705278 \n",
".. ... ... ... ... ... ... \n",
"145 156.579188 358.881977 5.045805 3.345408 2.849412 14.947433 \n",
"146 158.532617 221.374213 5.030398 3.244671 2.989694 17.112216 \n",
"147 160.672431 83.866449 4.976185 3.342438 2.972667 22.005111 \n",
"148 163.069275 306.358685 4.940651 3.287877 3.090450 20.656204 \n",
"149 180.000000 0.000000 4.731335 3.345559 3.345559 25.109194 \n",
" polar_ang azimuthal_ang phase_Vp phase_Vs1 phase_Vs2 group_Vp \\\n",
"0 0.0 0.0 4.731335 3.345559 3.345559 25.109194 \n",
"1 16.9 0.0 4.932846 3.345559 3.040613 18.736040 \n",
"2 19.3 222.5 4.988421 3.262362 3.040346 25.801391 \n",
"3 21.5 85.0 5.013496 3.343000 2.908654 21.006679 \n",
"4 23.4 307.5 5.067123 3.234312 2.938460 16.705278 \n",
".. ... ... ... ... ... ... \n",
"145 156.6 358.9 5.045805 3.345408 2.849412 14.947433 \n",
"146 158.5 221.4 5.030398 3.244671 2.989694 17.112216 \n",
"147 160.7 83.9 4.976185 3.342438 2.972667 22.005111 \n",
"148 163.1 306.4 4.940651 3.287877 3.090450 20.656204 \n",
"149 180.0 0.0 4.731335 3.345559 3.345559 25.109194 \n",
"\n",
" group_Vs1 group_Vs2 \n",
"0 17.754881 17.754881 \n",
Expand Down

0 comments on commit 0eb0b4f

Please sign in to comment.