Appendix 1: AT-Cut alpha-Quartz Resonator

Generated from QUARTZ.wxm with rfMaxima 0.2.5 on Maxima 5.46.0. All results were computed by executing the notebook; none is transcribed.

(See: Microsystem Design by Stephen D. Senturia, pp. 192-193, 570-577)
(See: An Introduction to the Theory of Piezoelectricity by Jiashi Yang, pp. 89-93)

Linear Piezoelectricity

Displacement (u), Strain (S), and Stress (T)

u: displacement [m]

u:matrix([u_x(x,y,z,t)],[u_y(x,y,z,t)],[u_z(x,y,z,t)])$
u2S(u):=matrix([diff(u[1,1],x)],[diff(u[2,1],y)],[diff(u[3,1],z)],[diff(u[2,1],z)+diff(u[3,1],y)],[diff(u[1,1],z)+diff(u[3,1],x)],[diff(u[1,1],y)+diff(u[2,1],x)])$

S: strain [] (%epsilon: normal strain, %gamma: shear strain)

S:matrix([%epsilon_x(x,y,z,t)],[%epsilon_y(x,y,z,t)],[%epsilon_z(x,y,z,t)],[%gamma_yz(x,y,z,t)],[%gamma_zx(x,y,z,t)],[%gamma_xy(x,y,z,t)])$
S=u2S(u);
(ϵx(x,y,z,t)ϵy(x,y,z,t)ϵz(x,y,z,t)γyz(x,y,z,t)γzx(x,y,z,t)γxy(x,y,z,t))=(dd x ux(x,y,z,t)dd y uy(x,y,z,t)dd z uz(x,y,z,t)dd y uz(x,y,z,t)+dd z uy(x,y,z,t)dd x uz(x,y,z,t)+dd z ux(x,y,z,t)dd x uy(x,y,z,t)+dd y ux(x,y,z,t))\begin{pmatrix} \epsilon_{x}\left(x , y , z , t\right)\\ \epsilon_{y}\left(x , y , z , t\right)\\ \epsilon_{z}\left(x , y , z , t\right)\\ \gamma_{yz}\left(x , y , z , t\right)\\ \gamma_{zx}\left(x , y , z , t\right)\\ \gamma_{xy}\left(x , y , z , t\right) \end{pmatrix} =\begin{pmatrix} {{d}\over{d\,x}}\,\mathit{u}_{x}\left(x , y , z , t\right)\\ {{d}\over{d\,y}}\,\mathit{u}_{y}\left(x , y , z , t\right)\\ {{d}\over{d\,z}}\,\mathit{u}_{z}\left(x , y , z , t\right)\\ {{d}\over{d\,y}}\,\mathit{u}_{z}\left(x , y , z , t\right)+{{d}\over{d\,z}}\,\mathit{u}_{y}\left(x , y , z , t\right)\\ {{d}\over{d\,x}}\,\mathit{u}_{z}\left(x , y , z , t\right)+{{d}\over{d\,z}}\,\mathit{u}_{x}\left(x , y , z , t\right)\\ {{d}\over{d\,x}}\,\mathit{u}_{y}\left(x , y , z , t\right)+{{d}\over{d\,y}}\,\mathit{u}_{x}\left(x , y , z , t\right) \end{pmatrix}
grad(u):=matrix([diff(u[1,1],x,1),diff(u[1,1],y,1),diff(u[1,1],z,1)],[diff(u[2,1],x,1),diff(u[2,1],y,1),diff(u[2,1],z,1)],[diff(u[3,1],x,1),diff(u[3,1],y,1),diff(u[3,1],z,1)])$
grad(u);
1/2*(grad(u)+transpose(grad(u)));
(dd x ux(x,y,z,t)dd y ux(x,y,z,t)dd z ux(x,y,z,t)dd x uy(x,y,z,t)dd y uy(x,y,z,t)dd z uy(x,y,z,t)dd x uz(x,y,z,t)dd y uz(x,y,z,t)dd z uz(x,y,z,t))\begin{pmatrix} {{d}\over{d\,x}}\,\mathit{u}_{x}\left(x , y , z , t\right)&{{d}\over{d\,y}}\,\mathit{u}_{x}\left(x , y , z , t\right)&{{d}\over{d\,z}}\,\mathit{u}_{x}\left(x , y , z , t\right)\\ {{d}\over{d\,x}}\,\mathit{u}_{y}\left(x , y , z , t\right)&{{d}\over{d\,y}}\,\mathit{u}_{y}\left(x , y , z , t\right)&{{d}\over{d\,z}}\,\mathit{u}_{y}\left(x , y , z , t\right)\\ {{d}\over{d\,x}}\,\mathit{u}_{z}\left(x , y , z , t\right)&{{d}\over{d\,y}}\,\mathit{u}_{z}\left(x , y , z , t\right)&{{d}\over{d\,z}}\,\mathit{u}_{z}\left(x , y , z , t\right) \end{pmatrix}(dd x ux(x,y,z,t)dd x uy(x,y,z,t)+dd y ux(x,y,z,t)2dd x uz(x,y,z,t)+dd z ux(x,y,z,t)2dd x uy(x,y,z,t)+dd y ux(x,y,z,t)2dd y uy(x,y,z,t)dd y uz(x,y,z,t)+dd z uy(x,y,z,t)2dd x uz(x,y,z,t)+dd z ux(x,y,z,t)2dd y uz(x,y,z,t)+dd z uy(x,y,z,t)2dd z uz(x,y,z,t))\begin{pmatrix} {{d}\over{d\,x}}\,\mathit{u}_{x}\left(x , y , z , t\right)&{{{{d}\over{d\,x}}\,\mathit{u}_{y}\left(x , y , z , t\right)+{{d}\over{d\,y}}\,\mathit{u}_{x}\left(x , y , z , t\right)}\over{2}}&{{{{d}\over{d\,x}}\,\mathit{u}_{z}\left(x , y , z , t\right)+{{d}\over{d\,z}}\,\mathit{u}_{x}\left(x , y , z , t\right)}\over{2}}\\ {{{{d}\over{d\,x}}\,\mathit{u}_{y}\left(x , y , z , t\right)+{{d}\over{d\,y}}\,\mathit{u}_{x}\left(x , y , z , t\right)}\over{2}}&{{d}\over{d\,y}}\,\mathit{u}_{y}\left(x , y , z , t\right)&{{{{d}\over{d\,y}}\,\mathit{u}_{z}\left(x , y , z , t\right)+{{d}\over{d\,z}}\,\mathit{u}_{y}\left(x , y , z , t\right)}\over{2}}\\ {{{{d}\over{d\,x}}\,\mathit{u}_{z}\left(x , y , z , t\right)+{{d}\over{d\,z}}\,\mathit{u}_{x}\left(x , y , z , t\right)}\over{2}}&{{{{d}\over{d\,y}}\,\mathit{u}_{z}\left(x , y , z , t\right)+{{d}\over{d\,z}}\,\mathit{u}_{y}\left(x , y , z , t\right)}\over{2}}&{{d}\over{d\,z}}\,\mathit{u}_{z}\left(x , y , z , t\right) \end{pmatrix}

T: stress, force per unit area [N/m^2] (%sigma: normall stress, %tau: shear stress)

T:matrix([%sigma_x(x,y,z,t)],[%sigma_y(x,y,z,t)],[%sigma_z(x,y,z,t)],[%tau_yz(x,y,z,t)],[%tau_zx(x,y,z,t)],[%tau_xy(x,y,z,t)]);
divT(T):=matrix([diff(T[1,1],x,1)+diff(T[6,1],y,1)+diff(T[5,1],z,1)],[diff(T[6,1],x,1)+diff(T[2,1],y,1)+diff(T[4,1],z,1)],[diff(T[5,1],x,1)+diff(T[4,1],y,1)+diff(T[3,1],z,1)])$
divT(T);
(σx(x,y,z,t)σy(x,y,z,t)σz(x,y,z,t)τyz(x,y,z,t)τzx(x,y,z,t)τxy(x,y,z,t))\begin{pmatrix} \sigma_{x}\left(x , y , z , t\right)\\ \sigma_{y}\left(x , y , z , t\right)\\ \sigma_{z}\left(x , y , z , t\right)\\ \tau_{yz}\left(x , y , z , t\right)\\ \tau_{zx}\left(x , y , z , t\right)\\ \tau_{xy}\left(x , y , z , t\right) \end{pmatrix}(dd z τzx(x,y,z,t)+dd y τxy(x,y,z,t)+dd x σx(x,y,z,t)dd z τyz(x,y,z,t)+dd x τxy(x,y,z,t)+dd y σy(x,y,z,t)dd x τzx(x,y,z,t)+dd y τyz(x,y,z,t)+dd z σz(x,y,z,t))\begin{pmatrix} {{d}\over{d\,z}}\,\tau_{zx}\left(x , y , z , t\right)+{{d}\over{d\,y}}\,\tau_{xy}\left(x , y , z , t\right)+{{d}\over{d\,x}}\,\sigma_{x}\left(x , y , z , t\right)\\ {{d}\over{d\,z}}\,\tau_{yz}\left(x , y , z , t\right)+{{d}\over{d\,x}}\,\tau_{xy}\left(x , y , z , t\right)+{{d}\over{d\,y}}\,\sigma_{y}\left(x , y , z , t\right)\\ {{d}\over{d\,x}}\,\tau_{zx}\left(x , y , z , t\right)+{{d}\over{d\,y}}\,\tau_{yz}\left(x , y , z , t\right)+{{d}\over{d\,z}}\,\sigma_{z}\left(x , y , z , t\right) \end{pmatrix}

Electric Displacement (D) and Electric Field (E)

D: electric displacement [C/m^2]

D:matrix([D_x(x,y,z,t)],[D_y(x,y,z,t)],[D_z(x,y,z,t)]);
/* Divergence of D is a SCALAR.  0.2.4 returned a 3x1 column of the three
   separate terms, which then made Gauss's law read as three equations. */
divD(D):=diff(D[1,1],x,1)+diff(D[2,1],y,1)+diff(D[3,1],z,1)$
divD(D);
(Dx(x,y,z,t)Dy(x,y,z,t)Dz(x,y,z,t))\begin{pmatrix} \mathit{D}_{x}\left(x , y , z , t\right)\\ \mathit{D}_{y}\left(x , y , z , t\right)\\ \mathit{D}_{z}\left(x , y , z , t\right) \end{pmatrix}dd z Dz(x,y,z,t)+dd y Dy(x,y,z,t)+dd x Dx(x,y,z,t){{d}\over{d\,z}}\,\mathit{D}_{z}\left(x , y , z , t\right)+{{d}\over{d\,y}}\,\mathit{D}_{y}\left(x , y , z , t\right)+{{d}\over{d\,x}}\,\mathit{D}_{x}\left(x , y , z , t\right)

E: electric field [V/m]

E:matrix([E_x(x,y,z,t)],[E_y(x,y,z,t)],[E_z(x,y,z,t)]);
(Ex(x,y,z,t)Ey(x,y,z,t)Ez(x,y,z,t))\begin{pmatrix} \mathit{E}_{x}\left(x , y , z , t\right)\\ \mathit{E}_{y}\left(x , y , z , t\right)\\ \mathit{E}_{z}\left(x , y , z , t\right) \end{pmatrix}

Linear Constitutive Tensors

The number of independent non-zero constants depends on the crystal symmetry. Y-cut (AT-, BT-, ST-cut) alpha-quartz is a material of crystal class 32.

%epsilon: dielectric permittivity tensor [C/(V.m)] (at constant strain)

declare([%epsilon_0,%epsilon_11,%epsilon_22,%epsilon_23,%epsilon_32,%epsilon_33],[constant,real,scalar])$
%epsilon:matrix([%epsilon_11,0,0],[0,%epsilon_22,%epsilon_23],[0,%epsilon_32,%epsilon_33]);
(ϵ11000ϵ22ϵ230ϵ32ϵ33)\begin{pmatrix} \epsilon_{11}&0&0\\ 0&\epsilon_{22}&\epsilon_{23}\\ 0&\epsilon_{32}&\epsilon_{33} \end{pmatrix}

C: elastic stiffness tensor [N/m^2] (at constant electric field)

declare([C_11,C_12,C_13,C_14,C_22,C_23,C_24,C_33,C_34,C_44,C_55,C_56,C_66],[constant,real,scalar])$
C:matrix([C_11, C_12, C_13, C_14, 0, 0],[C_12, C_22, C_23, C_24, 0, 0],[C_13, C_23, C_33, C_34, 0, 0],[C_14, C_24, C_34, C_44, 0, 0],[0, 0, 0, 0, C_55, C_56],[0, 0, 0, 0, C_56, C_66]);
(C11C12C13C1400C12C22C23C2400C13C23C33C3400C14C24C34C44000000C55C560000C56C66)\begin{pmatrix} C_{11}&C_{12}&C_{13}&C_{14}&0&0\\ C_{12}&C_{22}&C_{23}&C_{24}&0&0\\ C_{13}&C_{23}&C_{33}&C_{34}&0&0\\ C_{14}&C_{24}&C_{34}&C_{44}&0&0\\ 0&0&0&0&C_{55}&C_{56}\\ 0&0&0&0&C_{56}&C_{66} \end{pmatrix}

e: piezoelectricity tensor [C/m^2]

declare([e_11,e_12,e_13,e_14,e_25,e_26,e_35,e_36],[constant,real,scalar])$
e:matrix([e_11,e_12,e_13,e_14,0,0],[0,0,0,0,e_25,e_26],[0,0,0,0,e_35,e_36]);
(e11e12e13e14000000e25e260000e35e36)\begin{pmatrix} e_{11}&e_{12}&e_{13}&e_{14}&0&0\\ 0&0&0&0&e_{25}&e_{26}\\ 0&0&0&0&e_{35}&e_{36} \end{pmatrix}

Lippmann's Linear Constitutive Equations in Charge-Stress Form

T:C.S-transpose(e).E;
D:e.S+%epsilon.E;
(−e11 Ex(x,y,z,t)+C14 γyz(x,y,z,t)+C13 ϵz(x,y,z,t)+C12 ϵy(x,y,z,t)+C11 ϵx(x,y,z,t)−e12 Ex(x,y,z,t)+C24 γyz(x,y,z,t)+C23 ϵz(x,y,z,t)+C22 ϵy(x,y,z,t)+C12 ϵx(x,y,z,t)−e13 Ex(x,y,z,t)+C34 γyz(x,y,z,t)+C33 ϵz(x,y,z,t)+C23 ϵy(x,y,z,t)+C13 ϵx(x,y,z,t)−e14 Ex(x,y,z,t)+C44 γyz(x,y,z,t)+C34 ϵz(x,y,z,t)+C24 ϵy(x,y,z,t)+C14 ϵx(x,y,z,t)−e35 Ez(x,y,z,t)−e25 Ey(x,y,z,t)+C55 γzx(x,y,z,t)+C56 γxy(x,y,z,t)−e36 Ez(x,y,z,t)−e26 Ey(x,y,z,t)+C56 γzx(x,y,z,t)+C66 γxy(x,y,z,t))\begin{pmatrix} -e_{11}\,\mathit{E}_{x}\left(x , y , z , t\right)+C_{14}\,\gamma_{yz}\left(x , y , z , t\right)+C_{13}\,\epsilon_{z}\left(x , y , z , t\right)+C_{12}\,\epsilon_{y}\left(x , y , z , t\right)+C_{11}\,\epsilon_{x}\left(x , y , z , t\right)\\ -e_{12}\,\mathit{E}_{x}\left(x , y , z , t\right)+C_{24}\,\gamma_{yz}\left(x , y , z , t\right)+C_{23}\,\epsilon_{z}\left(x , y , z , t\right)+C_{22}\,\epsilon_{y}\left(x , y , z , t\right)+C_{12}\,\epsilon_{x}\left(x , y , z , t\right)\\ -e_{13}\,\mathit{E}_{x}\left(x , y , z , t\right)+C_{34}\,\gamma_{yz}\left(x , y , z , t\right)+C_{33}\,\epsilon_{z}\left(x , y , z , t\right)+C_{23}\,\epsilon_{y}\left(x , y , z , t\right)+C_{13}\,\epsilon_{x}\left(x , y , z , t\right)\\ -e_{14}\,\mathit{E}_{x}\left(x , y , z , t\right)+C_{44}\,\gamma_{yz}\left(x , y , z , t\right)+C_{34}\,\epsilon_{z}\left(x , y , z , t\right)+C_{24}\,\epsilon_{y}\left(x , y , z , t\right)+C_{14}\,\epsilon_{x}\left(x , y , z , t\right)\\ -e_{35}\,\mathit{E}_{z}\left(x , y , z , t\right)-e_{25}\,\mathit{E}_{y}\left(x , y , z , t\right)+C_{55}\,\gamma_{zx}\left(x , y , z , t\right)+C_{56}\,\gamma_{xy}\left(x , y , z , t\right)\\ -e_{36}\,\mathit{E}_{z}\left(x , y , z , t\right)-e_{26}\,\mathit{E}_{y}\left(x , y , z , t\right)+C_{56}\,\gamma_{zx}\left(x , y , z , t\right)+C_{66}\,\gamma_{xy}\left(x , y , z , t\right) \end{pmatrix}(ϵ11 Ex(x,y,z,t)+e14 γyz(x,y,z,t)+e13 ϵz(x,y,z,t)+e12 ϵy(x,y,z,t)+e11 ϵx(x,y,z,t)ϵ23 Ez(x,y,z,t)+ϵ22 Ey(x,y,z,t)+e25 γzx(x,y,z,t)+e26 γxy(x,y,z,t)ϵ33 Ez(x,y,z,t)+ϵ32 Ey(x,y,z,t)+e35 γzx(x,y,z,t)+e36 γxy(x,y,z,t))\begin{pmatrix} \epsilon_{11}\,\mathit{E}_{x}\left(x , y , z , t\right)+e_{14}\,\gamma_{yz}\left(x , y , z , t\right)+e_{13}\,\epsilon_{z}\left(x , y , z , t\right)+e_{12}\,\epsilon_{y}\left(x , y , z , t\right)+e_{11}\,\epsilon_{x}\left(x , y , z , t\right)\\ \epsilon_{23}\,\mathit{E}_{z}\left(x , y , z , t\right)+\epsilon_{22}\,\mathit{E}_{y}\left(x , y , z , t\right)+e_{25}\,\gamma_{zx}\left(x , y , z , t\right)+e_{26}\,\gamma_{xy}\left(x , y , z , t\right)\\ \epsilon_{33}\,\mathit{E}_{z}\left(x , y , z , t\right)+\epsilon_{32}\,\mathit{E}_{y}\left(x , y , z , t\right)+e_{35}\,\gamma_{zx}\left(x , y , z , t\right)+e_{36}\,\gamma_{xy}\left(x , y , z , t\right) \end{pmatrix}

Some of the most important resonator phenomena (e.g. acceleration sensitivity) are due to nonlinear effects. Quartz has numerous higher order constants, e.g., 14 third-order and 23 fourth-order elastic constants, as well as 16 third-order piezoelectric coefficients are known (Source: "Quartz Crystal Resonators and Oscillators for Frequency Control and Timing Applications - A Tutorial", by John R. Vig, November 2008.)

Boundary Value Problem

Assume a thickness-shear vibration of a AT-cut quartz plate around the Z-axis, excited by an electric field along the Y-axis. Assume the plate can not move along the Y- and Z-axis: u_y(y,t)=0, u_z(y,t)=0.

%epsilon_x(x,y,z,t):=0$
%epsilon_y(x,y,z,t):=0$
%epsilon_z(x,y,z,t):=0$
%gamma_yz(x,y,z,t):=0$
%gamma_zx(x,y,z,t):=0$  /* 0.2.4 wrote %gamma_xz, which is not the name used in S */
%gamma_xy(x,y,z,t):=diff(u_x(y,t),y,1)$

Assume the X- and Z-component of the electric field to be zero: E_x(x,y,z,t)=0, E_z(x,y,z,t)=0.

E:matrix([0],[E_y(y,t)],[0]);
T:ev(C.S-transpose(e).E);
D:ev(e.S+%epsilon.E);
(0Ey(y,t)0)\begin{pmatrix} 0\\ \mathit{E}_{y}\left(y , t\right)\\ 0 \end{pmatrix}(0000C56 (dd y ux(y,t))−e25 Ey(y,t)C66 (dd y ux(y,t))−e26 Ey(y,t))\begin{pmatrix} 0\\ 0\\ 0\\ 0\\ C_{56}\,\left({{d}\over{d\,y}}\,\mathit{u}_{x}\left(y , t\right)\right)-e_{25}\,\mathit{E}_{y}\left(y , t\right)\\ C_{66}\,\left({{d}\over{d\,y}}\,\mathit{u}_{x}\left(y , t\right)\right)-e_{26}\,\mathit{E}_{y}\left(y , t\right) \end{pmatrix}(0e26 (dd y ux(y,t))+ϵ22 Ey(y,t)e36 (dd y ux(y,t))+ϵ32 Ey(y,t))\begin{pmatrix} 0\\ e_{26}\,\left({{d}\over{d\,y}}\,\mathit{u}_{x}\left(y , t\right)\right)+\epsilon_{22}\,\mathit{E}_{y}\left(y , t\right)\\ e_{36}\,\left({{d}\over{d\,y}}\,\mathit{u}_{x}\left(y , t\right)\right)+\epsilon_{32}\,\mathit{E}_{y}\left(y , t\right) \end{pmatrix}

Potential Energy

Electrostatic Energy Density, w_e

w_e:expand(ev(1/2*transpose(E).D));
e26 Ey(y,t) (dd y ux(y,t))2+ϵ22 Ey(y,t)22{{e_{26}\,\mathit{E}_{y}\left(y , t\right)\,\left({{d}\over{d\,y}}\,\mathit{u}_{x}\left(y , t\right)\right)}\over{2}}+{{\epsilon_{22}\,\mathit{E}_{y}\left(y , t\right)^2}\over{2}}

Elastic Energy Density, w_m

w_m:expand(ev(1/2*transpose(S).C.S));
C66 (dd y ux(y,t))22{{C_{66}\,\left({{d}\over{d\,y}}\,\mathit{u}_{x}\left(y , t\right)\right)^2}\over{2}}

Total Energy Density

w_t:w_e+w_m;
C66 (dd y ux(y,t))22+e26 Ey(y,t) (dd y ux(y,t))2+ϵ22 Ey(y,t)22{{C_{66}\,\left({{d}\over{d\,y}}\,\mathit{u}_{x}\left(y , t\right)\right)^2}\over{2}}+{{e_{26}\,\mathit{E}_{y}\left(y , t\right)\,\left({{d}\over{d\,y}}\,\mathit{u}_{x}\left(y , t\right)\right)}\over{2}}+{{\epsilon_{22}\,\mathit{E}_{y}\left(y , t\right)^2}\over{2}}

Newton's Law of Motion for a Continuum [PDE in Quartz Volume]

Force density

ev(T);
f:ratsimp(divT(ev(T)));
(0000C56 (dd y ux(y,t))−e25 Ey(y,t)C66 (dd y ux(y,t))−e26 Ey(y,t))\begin{pmatrix} 0\\ 0\\ 0\\ 0\\ C_{56}\,\left({{d}\over{d\,y}}\,\mathit{u}_{x}\left(y , t\right)\right)-e_{25}\,\mathit{E}_{y}\left(y , t\right)\\ C_{66}\,\left({{d}\over{d\,y}}\,\mathit{u}_{x}\left(y , t\right)\right)-e_{26}\,\mathit{E}_{y}\left(y , t\right) \end{pmatrix}(C66 (d2d y2 ux(y,t))−e26 (dd y Ey(y,t))00)\begin{pmatrix} C_{66}\,\left({{d^2}\over{d\,y^2}}\,\mathit{u}_{x}\left(y , t\right)\right)-e_{26}\,\left({{d}\over{d\,y}}\,\mathit{E}_{y}\left(y , t\right)\right)\\ 0\\ 0 \end{pmatrix}

Viscosity and Q Factor

%eta: viscosity tensor [N/(m^2.s)] %tau_m is the motional time constant %omega_m is the mechanical resonant radian frequency Q_m is the mechanical Q factor

%tau_m:1/(%omega_m*Q_m);
%eta:%tau_m*C;
1ωm Qm{{1}\over{\omega_{m}\,\mathit{Q}_{m}}}(C11ωm QmC12ωm QmC13ωm QmC14ωm Qm00C12ωm QmC22ωm QmC23ωm QmC24ωm Qm00C13ωm QmC23ωm QmC33ωm QmC34ωm Qm00C14ωm QmC24ωm QmC34ωm QmC44ωm Qm000000C55ωm QmC56ωm Qm0000C56ωm QmC66ωm Qm)\begin{pmatrix} {{C_{11}}\over{\omega_{m}\,\mathit{Q}_{m}}}&{{C_{12}}\over{\omega_{m}\,\mathit{Q}_{m}}}&{{C_{13}}\over{\omega_{m}\,\mathit{Q}_{m}}}&{{C_{14}}\over{\omega_{m}\,\mathit{Q}_{m}}}&0&0\\ {{C_{12}}\over{\omega_{m}\,\mathit{Q}_{m}}}&{{C_{22}}\over{\omega_{m}\,\mathit{Q}_{m}}}&{{C_{23}}\over{\omega_{m}\,\mathit{Q}_{m}}}&{{C_{24}}\over{\omega_{m}\,\mathit{Q}_{m}}}&0&0\\ {{C_{13}}\over{\omega_{m}\,\mathit{Q}_{m}}}&{{C_{23}}\over{\omega_{m}\,\mathit{Q}_{m}}}&{{C_{33}}\over{\omega_{m}\,\mathit{Q}_{m}}}&{{C_{34}}\over{\omega_{m}\,\mathit{Q}_{m}}}&0&0\\ {{C_{14}}\over{\omega_{m}\,\mathit{Q}_{m}}}&{{C_{24}}\over{\omega_{m}\,\mathit{Q}_{m}}}&{{C_{34}}\over{\omega_{m}\,\mathit{Q}_{m}}}&{{C_{44}}\over{\omega_{m}\,\mathit{Q}_{m}}}&0&0\\ 0&0&0&0&{{C_{55}}\over{\omega_{m}\,\mathit{Q}_{m}}}&{{C_{56}}\over{\omega_{m}\,\mathit{Q}_{m}}}\\ 0&0&0&0&{{C_{56}}\over{\omega_{m}\,\mathit{Q}_{m}}}&{{C_{66}}\over{\omega_{m}\,\mathit{Q}_{m}}} \end{pmatrix}

Newton's Equation of Motion

%rho: mass density [kg/m^3]

f=%rho*diff(u,t,2);
(C66 (d2d y2 ux(y,t))−e26 (dd y Ey(y,t))00)=(ρ (d2d t2 ux(x,y,z,t))ρ (d2d t2 uy(x,y,z,t))ρ (d2d t2 uz(x,y,z,t)))\begin{pmatrix} C_{66}\,\left({{d^2}\over{d\,y^2}}\,\mathit{u}_{x}\left(y , t\right)\right)-e_{26}\,\left({{d}\over{d\,y}}\,\mathit{E}_{y}\left(y , t\right)\right)\\ 0\\ 0 \end{pmatrix} =\begin{pmatrix} \rho\,\left({{d^2}\over{d\,t^2}}\,\mathit{u}_{x}\left(x , y , z , t\right)\right)\\ \rho\,\left({{d^2}\over{d\,t^2}}\,\mathit{u}_{y}\left(x , y , z , t\right)\right)\\ \rho\,\left({{d^2}\over{d\,t^2}}\,\mathit{u}_{z}\left(x , y , z , t\right)\right) \end{pmatrix}
f_x:ev(f[1,1]);
f_x=%rho*diff(u_x(y,t),t,2);
C66 (d2d y2 ux(y,t))−e26 (dd y Ey(y,t))C_{66}\,\left({{d^2}\over{d\,y^2}}\,\mathit{u}_{x}\left(y , t\right)\right)-e_{26}\,\left({{d}\over{d\,y}}\,\mathit{E}_{y}\left(y , t\right)\right)C66 (d2d y2 ux(y,t))−e26 (dd y Ey(y,t))=ρ (d2d t2 ux(y,t))C_{66}\,\left({{d^2}\over{d\,y^2}}\,\mathit{u}_{x}\left(y , t\right)\right)-e_{26}\,\left({{d}\over{d\,y}}\,\mathit{E}_{y}\left(y , t\right)\right)=\rho\,\left({{d^2}\over{d\,t^2}}\,\mathit{u}_{x}\left(y , t\right)\right)

Gauss' Law [PDE in Quartz Volume]

divD(ev(D))=0;
e26 (d2d y2 ux(y,t))+ϵ22 (dd y Ey(y,t))=0e_{26}\,\left({{d^2}\over{d\,y^2}}\,\mathit{u}_{x}\left(y , t\right)\right)+\epsilon_{22}\,\left({{d}\over{d\,y}}\,\mathit{E}_{y}\left(y , t\right)\right)=0

Calculating 'diff(E_y(y,t),y,1)

-(e_26*('diff(u_x(y,t),y,2)))/%epsilon_22;
−e26 (d2d y2 ux(y,t))ϵ22-{{e_{26}\,\left({{d^2}\over{d\,y^2}}\,\mathit{u}_{x}\left(y , t\right)\right)}\over{\epsilon_{22}}}

Calculating E_y(y,t) B_1 is an integration constant.

declare([B_1],[constant,real,scalar])$
define(E_y(y,t),-(e_26*(diff(u_x(y,t),y,1)))/%epsilon_22-B_1);
Ey(y,t):=−e26 (dd y ux(y,t))ϵ22−B1\mathit{E}_{y}\left(y , t\right):=-{{e_{26}\,\left({{d}\over{d\,y}}\,\mathit{u}_{x}\left(y , t\right)\right)}\over{\epsilon_{22}}}-B_{1}

Boundary Conditions

Assume a solution u_x(y,t):=A_1*sin(%xi*y)*exp(%i*%omega_m*t) with %xi:sqrt(%rho/(C_66*(1+k_26^2)))*%omega
k_26: sqrt(e_26^2/(%epsilon_22*C_66))

declare([A_1],[constant,real,scalar])$
u_x(y):=A_1*sin(%xi*y);
ratsimp(diff(u_x(y),y,1));
ratsimp(diff(u_x(y),y,2));
ux(y):=A1 sin⁡(ξ y)\mathit{u}_{x}\left(y\right):=A_{1}\,\sin \left(\xi\,y\right)A1 ξ cos⁡(ξ y)A_{1}\,\xi\,\cos \left(\xi\,y\right)−A1 ξ2 sin⁡(ξ y)-A_{1}\,\xi^2\,\sin \left(\xi\,y\right)
define(T_xy(y),C_66*(1+k_26^2)*%th(2)+e_26*B_1);
D_y(y):=-%epsilon_22*B_1;
define(E_y(y),-(e_26*%th(4))/%epsilon_22-B_1);
define(f_x(y),C_66*(1+k_26^2)*%th(4));
/* Renamed from f_s_z/f_p_z: this is a thickness-SHEAR problem in x-y and the
   field is E_y.  0.2.4 differentiated E_z, which does not exist here. */
define(f_s_x(y),C_66*%th(5));
define(f_p_x(y),-e_26*(diff(E_y(y),y,1)));
Txy(y):=A1 C66 ξ (k262+1) cos⁡(ξ y)+B1 e26\mathit{T}_{xy}\left(y\right):=A_{1}\,C_{66}\,\xi\,\left(k_{26}^2+1\right)\,\cos \left(\xi\,y\right)+B_{1}\,e_{26}Dy(y):=(−ϵ22) B1\mathit{D}_{y}\left(y\right):=\left(-\epsilon_{22}\right)\,B_{1}Ey(y):=−A1 e26 ξ cos⁡(ξ y)ϵ22−B1\mathit{E}_{y}\left(y\right):=-{{A_{1}\,e_{26}\,\xi\,\cos \left(\xi\,y\right)}\over{\epsilon_{22}}}-B_{1}fx(y):=−A1 C66 ξ2 (k262+1) sin⁡(ξ y)\mathit{f}_{x}\left(y\right):=-A_{1}\,C_{66}\,\xi^2\,\left(k_{26}^2+1\right)\,\sin \left(\xi\,y\right)fs,x(y):=−A1 C66 ξ2 sin⁡(ξ y)\mathit{f}_{s,x}\left(y\right):=-A_{1}\,C_{66}\,\xi^2\,\sin \left(\xi\,y\right)fp,x(y):=−A1 e262 ξ2 sin⁡(ξ y)ϵ22\mathit{f}_{p,x}\left(y\right):=-{{A_{1}\,e_{26}^2\,\xi^2\,\sin \left(\xi\,y\right)}\over{\epsilon_{22}}}

Calculating the potential %Phi(y) B_2 is an integration constant.

declare([B_2],[constant,real,scalar])$
%Phi(y):=e_26*u_x(y)/%epsilon_22+B_1*y+B_2;
Φ(y):=e26 ux(y)ϵ22+B1 y+B2\Phi\left(y\right):={{e_{26}\,\mathit{u}_{x}\left(y\right)}\over{\epsilon_{22}}}+B_{1}\,y+B_{2}

%tau_xy(y=-g_0/2)=0, %tau_xy(y=g_0/2)=0, %Phi(y=g_0/2)-%Phi(y=-g_0/2)=V, with, g_0: thickness [m] V: applied voltage [V]

/* 0.2.4 solved this with V(t) on the right-hand side but then divided by V
   further down, so the two never cancelled. */
solve([ev(T_xy(g_0/2))=0,ev(%Phi(g_0/2)-%Phi(-g_0/2))=V],[A_1,B_1])$
A_1:ratsimp(ev(A_1,%));
B_1:ratsimp(ev(B_1,%th(2)));
ratsimp(ev(A_1,%xi=%pi/g_0));
ratsimp(ev(B_1,%xi=%pi/g_0));
−ϵ22 e26 Vϵ22 C66 ξ g0 cos⁡(ξ g02) k262−2 e262 sin⁡(ξ g02)+ϵ22 C66 ξ g0 cos⁡(ξ g02)-{{\epsilon_{22}\,e_{26}\,V}\over{\epsilon_{22}\,C_{66}\,\xi\,g_{0}\,\cos \left({{\xi\,g_{0}}\over{2}}\right)\,k_{26}^2-2\,e_{26}^2\,\sin \left({{\xi\,g_{0}}\over{2}}\right)+\epsilon_{22}\,C_{66}\,\xi\,g_{0}\,\cos \left({{\xi\,g_{0}}\over{2}}\right)}}ϵ22 C66 ξ V cos⁡(ξ g02) k262+ϵ22 C66 ξ V cos⁡(ξ g02)ϵ22 C66 ξ g0 cos⁡(ξ g02) k262−2 e262 sin⁡(ξ g02)+ϵ22 C66 ξ g0 cos⁡(ξ g02){{\epsilon_{22}\,C_{66}\,\xi\,V\,\cos \left({{\xi\,g_{0}}\over{2}}\right)\,k_{26}^2+\epsilon_{22}\,C_{66}\,\xi\,V\,\cos \left({{\xi\,g_{0}}\over{2}}\right)}\over{\epsilon_{22}\,C_{66}\,\xi\,g_{0}\,\cos \left({{\xi\,g_{0}}\over{2}}\right)\,k_{26}^2-2\,e_{26}^2\,\sin \left({{\xi\,g_{0}}\over{2}}\right)+\epsilon_{22}\,C_{66}\,\xi\,g_{0}\,\cos \left({{\xi\,g_{0}}\over{2}}\right)}}ϵ22 V2 e26{{\epsilon_{22}\,V}\over{2\,e_{26}}}00

Elastic Wave Equation and Christoffel Equation

(See: https://sepwww.stanford.edu/data/media/public/docs/sep92/hector1/paper_html/node2.html)

K:ev(1/k*matrix([k_x,0,0,0,k_z,k_y],[0,k_y,0,k_z,0,k_x],[0,0,k_z,k_y,k_x,0]),k_x=0,k_z=0);
k:ev(sqrt(k_x^2+k_y^2+k_z^2),k_x=0,k_z=0);
(00000kyk0kyk0000000kyk00)\begin{pmatrix} 0&0&0&0&0&{{\mathit{k}_{y}}\over{k}}\\ 0&{{\mathit{k}_{y}}\over{k}}&0&0&0&0\\ 0&0&0&{{\mathit{k}_{y}}\over{k}}&0&0 \end{pmatrix}∣ky∣\left| \mathit{k}_{y}\right|
M:ratsimp(k^2/%rho*K.subst(C_66*(1+k_26^2),C_66,C).transpose(K));
ratsimp(eigenvalues(M));
ratsimp(eigenvectors(M));
((C66 k262+C66) ky4ρ k2000C22 ky4ρ k2C24 ky4ρ k20C24 ky4ρ k2C44 ky4ρ k2)\begin{pmatrix} {{\left(C_{66}\,k_{26}^2+C_{66}\right)\,\mathit{k}_{y}^4}\over{\rho\,k^2}}&0&0\\ 0&{{C_{22}\,\mathit{k}_{y}^4}\over{\rho\,k^2}}&{{C_{24}\,\mathit{k}_{y}^4}\over{\rho\,k^2}}\\ 0&{{C_{24}\,\mathit{k}_{y}^4}\over{\rho\,k^2}}&{{C_{44}\,\mathit{k}_{y}^4}\over{\rho\,k^2}} \end{pmatrix}[[−C442−2 C22 C44+4 C242+C222 ky4+(−C44−C22) ky42 ρ k2,C442−2 C22 C44+4 C242+C222 ky4+(C44+C22) ky42 ρ k2,(C66 k262+C66) ky4ρ k2],[1,1,1]]\left[ \left[ -{{\sqrt{C_{44}^2-2\,C_{22}\,C_{44}+4\,C_{24}^2+C_{22}^2}\,\mathit{k}_{y}^4+\left(-C_{44}-C_{22}\right)\,\mathit{k}_{y}^4}\over{2\,\rho\,k^2}} , {{\sqrt{C_{44}^2-2\,C_{22}\,C_{44}+4\,C_{24}^2+C_{22}^2}\,\mathit{k}_{y}^4+\left(C_{44}+C_{22}\right)\,\mathit{k}_{y}^4}\over{2\,\rho\,k^2}} , {{\left(C_{66}\,k_{26}^2+C_{66}\right)\,\mathit{k}_{y}^4}\over{\rho\,k^2}} \right] , \left[ 1 , 1 , 1 \right] \right][[[−C442−2 C22 C44+4 C242+C222 ky4+(−C44−C22) ky42 ρ k2,C442−2 C22 C44+4 C242+C222 ky4+(C44+C22) ky42 ρ k2,(C66 k262+C66) ky4ρ k2],[1,1,1]],[[[0,1,−C442−2 C22 C44+4 C242+C222−C44+C222 C24]],[[0,1,C442−2 C22 C44+4 C242+C222+C44−C222 C24]],[[1,0,0]]]]\left[ \left[ \left[ -{{\sqrt{C_{44}^2-2\,C_{22}\,C_{44}+4\,C_{24}^2+C_{22}^2}\,\mathit{k}_{y}^4+\left(-C_{44}-C_{22}\right)\,\mathit{k}_{y}^4}\over{2\,\rho\,k^2}} , {{\sqrt{C_{44}^2-2\,C_{22}\,C_{44}+4\,C_{24}^2+C_{22}^2}\,\mathit{k}_{y}^4+\left(C_{44}+C_{22}\right)\,\mathit{k}_{y}^4}\over{2\,\rho\,k^2}} , {{\left(C_{66}\,k_{26}^2+C_{66}\right)\,\mathit{k}_{y}^4}\over{\rho\,k^2}} \right] , \left[ 1 , 1 , 1 \right] \right] , \left[ \left[ \left[ 0 , 1 , -{{\sqrt{C_{44}^2-2\,C_{22}\,C_{44}+4\,C_{24}^2+C_{22}^2}-C_{44}+C_{22}}\over{2\,C_{24}}} \right] \right] , \left[ \left[ 0 , 1 , {{\sqrt{C_{44}^2-2\,C_{22}\,C_{44}+4\,C_{24}^2+C_{22}^2}+C_{44}-C_{22}}\over{2\,C_{24}}} \right] \right] , \left[ \left[ 1 , 0 , 0 \right] \right] \right] \right]

The vibration frequency of the thickness-shear mode, f_m, can be used to determine the thickness of the Y-cut alpha-quartz plate, g_0

/* Without these, ratsimp and limit further down keep stopping to ask
   whether %epsilon_22 is positive, which makes the notebook impossible to
   batch and forces the reader to answer by hand.  0.2.4 declared
   only the four symbols on the first line. */
/* NB: %omega, %xi and %omega_m are deliberately NOT declared constant --
   limit(...,%omega,0) below needs %omega to be a variable.  The assume()
   calls are what silence the sign questions. */
declare([%rho,C_66,g_0,k_26,%epsilon_22,A,V],[constant,real,scalar])$
assume(%rho>0)$
assume(C_66>0)$
assume(g_0>0)$
assume(k_26>0)$
assume(%epsilon_22>0)$
assume(A>0)$
assume(%omega>0)$
assume(%xi>0)$
assume(%omega_m>0)$
assume(Q_m>0)$
k_y:2*%pi/(2*g_0);
%omega_m:sqrt(ev((k_y^2*C_66*(1+k_26^2))/(%rho)));
f_m:1/(2*%pi)*%omega_m;
πg0{{\pi}\over{g_{0}}}C66 π k262+1ρ g0{{\sqrt{C_{66}}\,\pi\,\sqrt{k_{26}^2+1}}\over{\sqrt{\rho}\,g_{0}}}C66 k262+12 ρ g0{{\sqrt{C_{66}}\,\sqrt{k_{26}^2+1}}\over{2\,\sqrt{\rho}\,g_{0}}}

m is the effective mass. A is the electrode area.

declare([A],[constant,real,scalar])$
/* Effective (modal) mass.  Two changes from 0.2.4:
     - integrate over the WHOLE plate, -g_0/2 .. g_0/2.  0.2.4 integrated
       over 0 .. g_0/2 only, i.e. half of it.
     - evaluate the mode shape at the fundamental thickness resonance,
       %xi = %pi/g_0, BEFORE integrating.  With %xi left symbolic Maxima
       stops and asks about floor(%xi*g_0/(4*%pi)) + ceiling(...) and never
       gets past it, which is one reason this notebook could not be batched.
   The answer is the familiar m = %rho*A*g_0/2, half the physical mass. */
m:A*%rho*integrate(ev(u_x(y)^2,%xi=%pi/g_0),y,-g_0/2,g_0/2)/A_1^2;
ρ A g0 (ϵ22 C66 ξ g0 cos⁡(ξ g02) k262−2 e262 sin⁡(ξ g02)+ϵ22 C66 ξ g0 cos⁡(ξ g02))22 (−2 e262 sin⁡(g0 ξ2)+ϵ22 C66 g0 k262 ξ cos⁡(g0 ξ2)+ϵ22 C66 g0 ξ cos⁡(g0 ξ2))2{{\rho\,A\,g_{0}\,\left(\epsilon_{22}\,C_{66}\,\xi\,g_{0}\,\cos \left({{\xi\,g_{0}}\over{2}}\right)\,k_{26}^2-2\,e_{26}^2\,\sin \left({{\xi\,g_{0}}\over{2}}\right)+\epsilon_{22}\,C_{66}\,\xi\,g_{0}\,\cos \left({{\xi\,g_{0}}\over{2}}\right)\right)^2}\over{2\,\left(-2\,e_{26}^2\,\sin \left({{g_{0}\,\xi}\over{2}}\right)+\epsilon_{22}\,C_{66}\,g_{0}\,k_{26}^2\,\xi\,\cos \left({{g_{0}\,\xi}\over{2}}\right)+\epsilon_{22}\,C_{66}\,g_{0}\,\xi\,\cos \left({{g_{0}\,\xi}\over{2}}\right)\right)^2}}

k_1_eff is the linear spring constant.

solve([%omega_m=sqrt(k_1_eff/m)],[k_1_eff])$
/* 0.2.4 evaluated this but never assigned it, so %eta below used an unbound
   k_1_eff. */
k_1_eff:ev(k_1_eff,%);
π2 A C66 k262+π2 A C662 g0{{\pi^2\,A\,C_{66}\,k_{26}^2+\pi^2\,A\,C_{66}}\over{2\,g_{0}}}

Tf is the linear temperature coefficient of frequency.

/* Placeholder: nothing in f_m depends on Temp yet, so this evaluates to 0.
   Give the stiffness, density and thickness a temperature dependence to get
   a meaningful turnover curve. */
Tf:'diff(log(f_m),Temp);
dd Temp log⁡(C66 k262+12 ρ g0){{d}\over{d\,\mathit{Temp}}}\,\log \left({{\sqrt{C_{66}}\,\sqrt{k_{26}^2+1}}\over{2\,\sqrt{\rho}\,g_{0}}}\right)

Charge, q, on the Electrodes

Boundary condition for D A is the electrode area.

q:ratsimp(-D_y(g_0/2)*A);
−(ϵ222 A C66 V k262+ϵ222 A C66 V) ξ cos⁡(g0 ξ2)2 e262 sin⁡(g0 ξ2)+(−ϵ22 C66 g0 k262−ϵ22 C66 g0) ξ cos⁡(g0 ξ2)-{{\left(\epsilon_{22}^2\,A\,C_{66}\,V\,k_{26}^2+\epsilon_{22}^2\,A\,C_{66}\,V\right)\,\xi\,\cos \left({{g_{0}\,\xi}\over{2}}\right)}\over{2\,e_{26}^2\,\sin \left({{g_{0}\,\xi}\over{2}}\right)+\left(-\epsilon_{22}\,C_{66}\,g_{0}\,k_{26}^2-\epsilon_{22}\,C_{66}\,g_{0}\right)\,\xi\,\cos \left({{g_{0}\,\xi}\over{2}}\right)}}

Motional Capacitance, C_m, Static Capacitance, C_0, and Transduction Factor, %eta

C_XTAL:ev(q/V,
    %xi=sqrt(%rho/(C_66*(1+e_26^2/(%epsilon_22*C_66))))*%omega,
    k_26=sqrt(e_26^2/(%epsilon_22*C_66)));
C_0:ratsimp(limit(C_XTAL,%omega,0));
ratsimp(C_XTAL-C_0)$
C_m:ratsimp(ev(k_26^2*C_0,
    k_26=sqrt(e_26^2/(%epsilon_22*C_66))));
/* Transduction factor: %eta = sqrt(C_m * k_1_eff), with C_m = k_26^2 C_0.
   FBAR.wxm carried an extra (1+k^2) factor; the two notebooks now agree. */
%eta:ratsimp(ev(sqrt(C_m*k_1_eff),
    k_26=sqrt(e_26^2/(%epsilon_22*C_66))));
−ρ (ϵ22 A V e262+ϵ222 A C66 V) ω cos⁡(ρ g0 ω2 C66 e262ϵ22 C66+1)C66 V e262ϵ22 C66+1 (2 e262 sin⁡(ρ g0 ω2 C66 e262ϵ22 C66+1)+ρ (−e262 g0−ϵ22 C66 g0) ω cos⁡(ρ g0 ω2 C66 e262ϵ22 C66+1)C66 e262ϵ22 C66+1)-{{\sqrt{\rho}\,\left(\epsilon_{22}\,A\,V\,e_{26}^2+\epsilon_{22}^2\,A\,C_{66}\,V\right)\,\omega\,\cos \left({{\sqrt{\rho}\,g_{0}\,\omega}\over{2\,\sqrt{C_{66}}\,\sqrt{{{e_{26}^2}\over{\epsilon_{22}\,C_{66}}}+1}}}\right)}\over{\sqrt{C_{66}}\,V\,\sqrt{{{e_{26}^2}\over{\epsilon_{22}\,C_{66}}}+1}\,\left(2\,e_{26}^2\,\sin \left({{\sqrt{\rho}\,g_{0}\,\omega}\over{2\,\sqrt{C_{66}}\,\sqrt{{{e_{26}^2}\over{\epsilon_{22}\,C_{66}}}+1}}}\right)+{{\sqrt{\rho}\,\left(-e_{26}^2\,g_{0}-\epsilon_{22}\,C_{66}\,g_{0}\right)\,\omega\,\cos \left({{\sqrt{\rho}\,g_{0}\,\omega}\over{2\,\sqrt{C_{66}}\,\sqrt{{{e_{26}^2}\over{\epsilon_{22}\,C_{66}}}+1}}}\right)}\over{\sqrt{C_{66}}\,\sqrt{{{e_{26}^2}\over{\epsilon_{22}\,C_{66}}}+1}}}\right)}}A e262+ϵ22 A C66C66 g0{{A\,e_{26}^2+\epsilon_{22}\,A\,C_{66}}\over{C_{66}\,g_{0}}}A e264+ϵ22 A C66 e262ϵ22 C662 g0{{A\,e_{26}^4+\epsilon_{22}\,A\,C_{66}\,e_{26}^2}\over{\epsilon_{22}\,C_{66}^2\,g_{0}}}π2 A e262+ϵ22 π2 A C66 A e264+ϵ22 A C66 e2622 ϵ22 C66 g0{{\sqrt{\pi^2\,A\,e_{26}^2+\epsilon_{22}\,\pi^2\,A\,C_{66}}\,\sqrt{A\,e_{26}^4+\epsilon_{22}\,A\,C_{66}\,e_{26}^2}}\over{\sqrt{2}\,\epsilon_{22}\,C_{66}\,g_{0}}}

Material Parameters

AT-Cut alpha-Quartz (Trigonal Structure, Crystal Class 32)

(See "Elastic and Piezoelectric Constants of Alpha-Quartz", by R. Bechmann, American Physical Society Phys. Rev., vol. 110, no. 5, pp. 1060-1061, June 1958)

C_11:86.74;
C_12:-8.25;
C_13:27.15;
C_14:-3.66;
C_22:129.77;
C_23:-7.42;
C_24:5.70;
C_33:102.83;
C_34:9.92;
C_44:38.61;
C_55:68.81;
C_56:2.53;
C_66:29.01;
C:1e9*matrix([C_11, C_12, C_13, C_14, 0, 0],[C_12, C_22, C_23, C_24, 0, 0],[C_13, C_23, C_33, C_34, 0, 0],[C_14, C_24, C_34, C_44, 0, 0],[0, 0, 0, 0, C_55, C_56],[0, 0, 0, 0, C_56, C_66]);
86.7486.74−8.25-8.2527.1527.15−3.66-3.66129.77129.77−7.42-7.425.75.7102.83102.839.929.9238.6138.6168.8168.812.532.5329.0129.01(8.674×10+10−8.25×10+92.715×10+10−3.66×10+900−8.25×10+91.2977×10+11−7.42×10+95.7×10+9002.715×10+10−7.42×10+91.0283×10+119.92×10+900−3.66×10+95.7×10+99.92×10+93.861×10+100000006.881×10+102.53×10+900002.53×10+92.901×10+10)\begin{pmatrix} 8.674 \times 10^{+10}&-8.25 \times 10^{+9}&2.715 \times 10^{+10}&-3.66 \times 10^{+9}&0&0\\ -8.25 \times 10^{+9}&1.2977 \times 10^{+11}&-7.42 \times 10^{+9}&5.7 \times 10^{+9}&0&0\\ 2.715 \times 10^{+10}&-7.42 \times 10^{+9}&1.0283 \times 10^{+11}&9.92 \times 10^{+9}&0&0\\ -3.66 \times 10^{+9}&5.7 \times 10^{+9}&9.92 \times 10^{+9}&3.861 \times 10^{+10}&0&0\\ 0&0&0&0&6.881 \times 10^{+10}&2.53 \times 10^{+9}\\ 0&0&0&0&2.53 \times 10^{+9}&2.901 \times 10^{+10} \end{pmatrix}
e_11:0.171;
e_12:-0.152;
e_13:-0.0187;
e_14:0.067;
e_25:0.108;
e_26:-0.095;
e_35:-0.0761;
e_36:0.067;
e:matrix([e_11,e_12,e_13,e_14,0,0],[0,0,0,0,e_25,e_26],[0,0,0,0,e_35,e_36]);
0.1710.171−0.152-0.152−0.0187-0.01870.0670.0670.1080.108−0.095-0.095−0.0761-0.07610.0670.067(0.171−0.152−0.01870.0670000000.108−0.0950000−0.07610.067)\begin{pmatrix} 0.171&-0.152&-0.0187&0.067&0&0\\ 0&0&0&0&0.108&-0.095\\ 0&0&0&0&-0.0761&0.067 \end{pmatrix}
%epsilon_11:39.21;
%epsilon_22:39.82;
%epsilon_23:0.86;
%epsilon_32:0.86;
%epsilon_33:40.42;
%epsilon:1e-12*matrix([%epsilon_11,0,0],[0,%epsilon_22,%epsilon_23],[0,%epsilon_32,%epsilon_33]);
39.2139.2139.8239.820.860.860.860.8640.4240.42(3.921×10−110003.982×10−118.6×10−1308.6×10−134.042×10−11)\begin{pmatrix} 3.921 \times 10^{-11}&0&0\\ 0&3.982 \times 10^{-11}&8.6 \times 10^{-13}\\ 0&8.6 \times 10^{-13}&4.042 \times 10^{-11} \end{pmatrix}
%rho:2650.67;
2650.672650.67

k: piezoelectric coupling coefficient

k_26:sqrt(e_26^2/(C_66*1e9*%epsilon_22*1e-12));
k_26^2;
0.088389181656823840.088389181656823840.0078126474339630050.007812647433963005