Message Boards Message Boards

0
|
3486 Views
|
4 Replies
|
1 Total Likes
View groups...
Share
Share this post:

Graphical Representation of the Energy Surface of a Triaxial Rotor

Energy surface - Mathematica

I am working with a triaxial rigid rotor, which is a 3D object with no symmetry axis, so the moments of inertia about the principal axes are all different. It is well known from textbooks that the energy (rotational) of a triaxial rotor is given by the formula: $$E=A_1I_1^2+A_2I_2^2+A_3I_3^2\ ,$$ where the terms $A_k=\frac{1}{2\mathcal{I}_k}$ are the inverse of the moments of inertia $\mathcal{I}_k$, and $I_k$ represent the components of the total angular momentum $\mathbf{I}$.
My problem involves plotting the energy surface of the triaxial rotor, which should be a triaxial ellipsoid, deformed accordingly to the values of the inertia parameters $A_k$.

I chose to work in spherical coordinates, so a representation of the angular momentum vector using its components will be something like this: $$I_1\equiv x_1 =r\sin\theta\cos\varphi \\ I_2\equiv x_2=r\sin\theta\sin\varphi \\ I_3\equiv x_3=r\cos\theta$$ where $r\in[0,I]$ (can be any number but no larger than the maximum spin value $I$) , $\theta\in[0,\pi]$ (polar angle) and $\varphi\in[0,2\pi]$ (azimuthal angle).

The following arbitrary values for the inertia factors $A3=1$ , $A_1=6$, $A_2=3$, and the maximum spin $I=2$ are given. Note that the smallest inertia factor is the one corresponding to the 3-axis (so it means that the 3-axis has the largest moment of inertia). Therefore, choosing $x_3$ as the quantization axis is the right choice.

Now, for the graphical representation of the energy ellipsoid $E$, I was thinking of doing a SphericalPlot3D, since I work in spherical coordinates. Fixing $r$ to a value of $I=2$, and then trying to plot $E=f(r=I,\theta,\varphi)$, should (?) give me an ellipsoid. However, this is not what I obtain. The 3d object has some kind of neck. I don't understand why that is.

This is the code:

x1[r_, theta_, fi_] := r Sin[theta] Cos[fi]; 
x2[r_, theta_, fi_] := r Sin[theta] Sin[fi]; 
x3[r_, theta_, fi_] := r Cos[theta];

EnSphere[a1_, a2_, a3_, r_, theta_, fi_] := 
  a1*x1[r, theta, fi]^2 + a2*x2[r, theta, fi]^2 + 
   a3*x3[r, theta, fi]^2;

SphericalPlot3D[
 En3[A1, A2, A3, 1, \[Theta], \[Phi]], {\[Theta], 
  0, \[Pi]}, {\[Phi], 0, 2\[Pi]}, PlotStyle -> Blue]

This is what I obtain

image1

This is what I expected to obtain

image2

Why is this happening? And how can I solve this? Is there a way to obtain the 3D energy surface of the triaxial rigid rotor starting from the values of the inertial parameters $A_k$?

Thank you in advance!

POSTED BY: Robert Poenaru
4 Replies

Well, I think we have two totally different representations. The contourplot asks for all momenta {x,y,z} with a given energie (1.6 here , Gianluca looked at energy = 1 )

ContourPlot3D[ 1.6 == 6 x^2 + 3 y^2 + z^2, {x, -1, 1}, {y, -1, 1}, {z, -1.3, 1.3}]

(or

RegionPlot3D[ 0 < 6 x^2 + 3 y^2 + z^2 < 1.6, {x, -1, 1}, {y, -1, 1}, {z, -1, 1}]

)

But the momenta giving this energy don't have all the length 1. E.g. chose a momentum { .2, .3, z } and ask, for which z the energy 1.6 is achieved . You will get a solution, but this has not the length 1.

Clear[x, y, z]
vec = {.2, .3, z} /. Solve[1.6 == 6 x^2 + 3 y^2 + z^2 /. x -> .2 /. y -> .3, z][[1, 1]]
vec.vec

If on the other hand you ask for the energies produced by a momentum of length 1 you have that these energies are angle-dependent.

Plot[6 x1[1, Pi/2, \[Phi]]^2 + 3 x2[1, Pi/2, \[Phi]]^2 +   1 x3[1, Pi/2, \[Phi]]^2, {\[Phi], 0, 2 \[Pi]}]
Plot[6 x1[1, th, 0]^2 + 3 x2[1, th, 0]^2 + 1 x3[1, th, 0]^2, {th, 0,  Pi}]

This explains the first figure in blue: it is the energy given by momenta of length 1 pointing in different directions.

POSTED BY: Hans Dolhaine

SphericalPlot3D works as a parametric plot, not as an implicit plot. You must first solve with respect to r and then feed the result to SphericalPlot3D:

Block[{a1 = 1, a2 = 2, a3 = 3},
 sol = Solve[EnSphere[a1, a2, a3, r, \[Theta], \[Phi]] == 1, r];
 rAsFunctionOFThetaFi = r /. sol[[2]];
 SphericalPlot3D[
  Evaluate[rAsFunctionOFThetaFi], {\[Theta], 0, \[Pi]}, {\[Phi], 0, 
   2 \[Pi]}, PlotStyle -> Blue]]
POSTED BY: Gianluca Gorni

Like this perhaps?

ContourPlot3D[ a1 x^2 + a2 y^2 + a3 z^2 /. {a1 -> 6, a2 -> 3, a3 -> 1}, 
{x, -1, 1}, {y, -1, 1}, {z, -1.3, 1.3}, Contours -> 1.5]
POSTED BY: Hans Dolhaine

Well, you used cartesian coordinates. If I switch to the components x1, x2, and x3 which are expressed in terms of spherical coordinates, then the result I get is different.

ContourPlot3D[
 a1 x1[r, \[Theta], \[Phi]]^2 + a2 x2[r, \[Theta], \[Phi]]^2 + 
   a3 x3[r, \[Theta], \[Phi]]^2 /. {a1 -> 6, a2 -> 3, a3 -> 1}, {r, 0,
   2}, {\[Theta], 0, \[Pi]}, {\[Phi], 0, 2 \[Pi]}, Contours -> 1.5]
POSTED BY: Robert Poenaru
Reply to this discussion
Community posts can be styled and formatted using the Markdown syntax.
Reply Preview
Attachments
Remove
or Discard

Group Abstract Group Abstract