NanoComp / meep

free finite-difference time-domain (FDTD) software for electromagnetic simulations
GNU General Public License v2.0
1.2k stars 612 forks source link

Apparent problem with geometry_lattice in MPB python tutorial #757

Closed JohnWeiner closed 5 years ago

JohnWeiner commented 5 years ago

In the section, Bands of a Triangular Lattice, in the MPB Python tutorial, the specification for the k_points is ms.k_points = [mp.Vector3(), # Gamma mp.Vector3(y=0.5), # M mp.Vector3(-1 / 3, 1 / 3), # K mp.Vector3()] # Gamma The K k_point, (-⅓, ⅓) does not correspond to the K point in the first Brillouin Zone (or the irreducible Brillouin Zone) as shown in Fig. 3, p. 237 or inset of Fig. 10, p. 76 of Photonic Crystals, Molding the Flow of Light, 2nd ed.. It is difficult to understand why x,y = -⅓,⅓ are specified.

JohnWeiner commented 5 years ago

The coefficients must correspond to the expansion coefficients for the k_points using the basis vectors of the triangular rod periodic structure. In the tutorial the angle between the basis vectors is 60 degrees. In the book, Photonic Crystals... the basis vectors are separated by 120 degrees (cf Fig. 3a, p. 237). I'm not sure yet if this difference is really important or not.

soamaven commented 5 years ago

The points in the book seem to be cartesian. If you do the conversion, you'll see they are the same in the MPB lattice basis of the tutorial.

JohnWeiner commented 5 years ago

I'm reopening this issue because I still don't get it. The mpb_tri_rods.py tutorial posits two basis vectors in reciprocal space (not real space). They are b1 = sqrt(3)/2 x + ½ y and b2 = sqrt(3)/2 x - ½ y. The **x,y are cartesian basis vectors. One can specify the coordinates of the Γ,M,K,Γ points of the irreducible Brillouin zone using either the cartesian basis or the b1,b2 basis. Looking at the left inset of Fig. 10, Chapter 5 (p. 76) of Photonic Crystals...2nd ed., the coordinates of the three points of the Brillouin zone appear to be (in **x,y basis): Γ = 0,0 M = 0,1/2 K = 1/(2sqrt(3)), ½ In order to use the two basis vectors b1,b2 to express the Γ,M,K points, it is convenient to rotate them ccw by 120 deg. Then b2 is oriented along the +y axis. The two rotated b1,b2 basis vectors, expanded in the cartesian basis, are then: b1 = -sqrt(3)/2 x + ½ y and b2= 0 x + 1 y. Now the coefficients for the Γ, M, K points in the rotated (b1,b2) basis should be: Γ = 0,0 M = 0,1/2 K = -⅓,+⅔ However, the tutorial specifies K = -⅓, +⅓ for the two coefficients, presumably in the (b1,b2) basis. It appears to me that a ⅓ coefficient for b2 will not give the correct y coordinate for K. The correct y value should be y = ½.
Any help with this issue would be greatly appreciated.

soamaven commented 5 years ago

Edit: Just see @stevengj 's answer below and ignore this one!

I wasn't clear before: if you use MPB/MEEP's Vector3.reciprocal_to_cartesian() method, they will have the same magnitude. I guess I never thought hard about what was going on before. But I still think it is correct.

I started to draw this out and I think I see your point.

I think the key thing is that we are trying to force a cartesian irreproducible brillouin zone directly overlayed onto the reciprocal triangular lattice. Something to keep in mind is that the Γ,M,K,Γ points in the book are indeed given in cartesian coordinates, and they are the irreproducible zone of a supercell, which is not the same thing we are dealing with. So trying to get them to match up is not happening because the change from real to reciprocal messes it up. Hopefully this diagram helps. I think I messed up the x and y but its too late now, and I am 90% sure that I am making the same mistake of overlaying real and reciprocal coordinates. The change of space when going from a supercell makes things harder. I did the best with my artistic abilities to draw out how, conceptionally, the MPB example is correct because the ultimate goal is to only traverse the irreproducible symmetry zone, of which there are 16 in triangular lattice rather than the 8 of the square lattice. You can see how the irreproducible zone in the supercell is just a scaled up version of the irreproducible zones on the right hand side, so the points are equivalent, but the zone is degenerate.

hexagonal

JohnWeiner commented 5 years ago

I was confused about the basis vectors in real space and in reciprocal space for quite some time, but now I think I am clear about which vectors belong to which space. For any discussion of k-points, Γ M K, the basis vectors must be in reciprocal space. Please take a look at Fig. 3, p. 237, Appendix B of Photonic Crystals...2nd ed.. The panel on the left shows basis vectors in real space. The panel in the middle shows the corresponding basis vectors in reciprocal space. The panel on the right shows the first Brillouin zone (shaded in yellow) in reciprocal space. The irreducible Brillouin zone is a sector as shown (shaded in blue) in the left inset of Fig.10, p. 76, Chapter 5 of the book. There are 12 equivalent irreducible Brillouin zones within the first Brillouin zone of the triangular lattice.
The two basis vectors, (b1,b2) in the center panel of Fig. 3 p. 237) are in reciprocal space and are oriented as in the MBP tutorial with a 60 degree separation. However these two basis vectors are not identical to (b1,b2) in the tutorial. The basis vectors in the book are (dropping the common factor of 2π/a): b1 = x + 1/sqrt(3) y and b2 = x - 1/sqrt(3) y. The basis vectors in the tutorial are: b1 = sqrt(3)/2 x + ½ y and b2 = sqrt(3)/2 x - ½ y. The difference is in the choice of length of the basis vectors. The lengths of (b1,b2) in the book are 2/sqrt(3). The lengths of (b1,b2) in the tutorial are 1. The difference in length is simply a choice of normalisation. However, I still don't understand the determination of (-⅓,⅓) for the coefficients of (b1,b2) that specify the k-point K in the tutorial. The coefficients (0,0), (0,1/2) for points Γ,M respectively do make sense. But it seems to me that one needs (-⅓,⅔) for a correct specification of K, whose coordinates are (in the Cartesian basis x,y of reciprocal space), (1/(2sqrt(3)),½). I still don't understand where (-⅓,⅓) comes from, and any enlightenment would be greatly appreciated.

stevengj commented 5 years ago

I haven't checked this recently, but I recall the irreducible Brillouin zone in the tutorial is a rotation of the one in the appendix of Photonic Crystals. (There are 12 possible choices of the irreducible Brillouin zone that are equivalent under symmetry operations for C₆ᵥ symmetry.)

stevengj commented 5 years ago

Something to keep in mind is that the Γ,M,K,Γ points in the book are indeed given in cartesian coordinates, and they are the irreproducible zone of a supercell, which is not the same thing we are dealing with

No, there is no supercell in the book: the given basis vectors (a₁,a₂ in real space or b₁,b₂ in reciprocal spacce) are primitive lattice vectors — the edges of a primitive cell (a rhombus): image Notice, however, that the appendix of the book uses a different choice of primitive lattice vectors a₁,a₂ than in the MPB triangular-lattice tutorial — in the tutorial, the whole lattice is rotated by 30° compared to the appendix, and the basis vectors (√3/2, ±1/2) have an angle of 60° between them (vs. 120° for the a₁,a₂ vectors in the appendix).

Realize that in any triangular lattice there are multiple valid choices of primitive lattice vectors (unit cells) and multiple valid definitions of the irreducible Brillouin zone.

soamaven commented 5 years ago

Thanks for the correction!

stevengj commented 5 years ago

Closing this as this seems to be more of an question of understanding the different coordinate systems and basis choices than an actual problem with the tutorial.

JohnWeiner commented 5 years ago

No, it is a real problem with the tutorial because the coefficients -⅓, ⅓ for the K k-point, which are simply posited in the tutorial, cannot be verified based on the (b1,b2) basis vector definition in the tutorial or in the book. I show in the earlier comment that the 120 deg rotation is necessary to make any sense of the coefficients. Please read carefully my comments starting from the one where I reopened the issue a day or so ago It is true that the basis vectors in reciprocal space (b1,b2) must be rotated by 120 deg (ccw) before the coefficients (0,0), (0,1/2) for k-points Γ,M respectively make sense. But consistent with this rotation are the coefficients (-⅓,⅔) NOT (-⅓,⅓) for the k-point K. The reasoning is spelled out in my two comments above. There is something about the choice of (-⅓ ,⅓) that makes it "correct" but whatever it is, it does not appear to come from the definition of the basis vectors in real space or reciprocal space. I agree with your comments about the basis vectors in real space and in reciprocal space and their choice, comparing the Fig in the book and the choice of b1,b2 in the tutorial. I tried to make this very clear in my comments above which are consistent with yours. I think it is important to get to the bottom of this issue and not simply accept coefficients listed in the tutorial as correct without verification. Sorry to be so insistent, and maybe I have just overlooked something that is obvious to everybody else. If so, I would like to know what it is.

soamaven commented 5 years ago

Let me know if this should be discussed elsewhere as you've already closed this.

No, there is no supercell in the book: the given basis vectors (a₁,a₂ in real space or b₁,b₂ in reciprocal spacce) are primitive lattice vectors — the edges of a primitive cell (a rhombus): image

@stevengj as a point of information, is not the image of the 6 circles in a triangular lattice representative of a super cell (see the first drawing of mine to the left and also pg 76, fig10 inset)? I think it is drawn this way because it is more intuitive/demonstrative from an educational stand-point, whereas the true primitive triangular cell would only include a single circle, and has the rhombus shape. So why does the Kappa point extend from one origin to the next if it is not a supercell? Should it not be only half-way between the two, otherwise it is redundant by reflection (see included image)? It is my understanding that this is the source of band-folding when calculating band-structures with FDTD, per Section 4.6 Adv. In Comp. FDTD. In the figure you posted, the far-right plot seems to support this; it shows the edge of the Brillouin zone at the bisect of the lines connecting the origins, as is correct. Just to be clear, my issue is not stemming from the image you posted from the Appendix, I'm on board with that, but the one in fig10 inset. Am i thinking in the wrong space-coordinate system?

SuperPrimitive

I'll admit to being wrong about the number of irriducible cells and my drawing of the "brillouin zone", I should probably not be responding so early into the mornings :/ Perhaps I have gotten so used to thinking in context of required supercell of the Cartesian coordinates in FDTD, that everything is a supercell :joy: .

I think the issue @JohnWeiner is having comes from the fact that the book and tutorial are meant to be equivalent in result, but were not meant to be equivalent in method? Just as there are multiple ways to approach a mathematical problem that will give the same correct result.

oskooi commented 5 years ago

For reference, here is a post from @stevengj to the mpb-discuss list from 8/27/2002 with subject Re: Triangular Rods Example (unfortunately, this post is not available from the archive which only goes back to 03/29/2006 but is accessible from a newsreader via the NNTP interface address: news.gmane.org/gmane.comp.science.photonic-bands):

On Tue, 27 Aug 2002, Todd Tolliver wrote:
> I have a simple question about the triangular rods example in the
> tutorial section of the MIT Photonic Bands Manual. In that example, on
> page 19, the vector describing the k-point for symmetry point 'K' is
> specified as (/ -3, / 3, 0). When looking at the irreducible Brillion
> zone, it seems that it is not possible for this vector to point in the
> correct direction. In other words, the vector be something like (/ -3,
> / 2, 0) instead. I realize that the definition of the k-point is in
> the basis of the reciprocal lattice vectors. Still a vector with equal
> x and y components would give the incrrect angle and therefore point
> to the wrong k-point. Any clarification would be greatly appreciated.

K = (-1/3, 1/3) is correct.  You are perhaps confused for three reasons.  
First, in the triangular lattice both the lattice and reciprocal bases are
non-orthogonal (affine), which makes them a little less intuitive than you
might like.  Second, there are many possible pairs of reciprocal lattice
basis vectors, and MPB's choice (fixed by the orthogonality relation with
the real-space lattice) may not be the same as the one you are drawing on
your sheet of hex paper.  Third, the reciprocal lattice is rotated by 30
degrees from the real-space lattice (see e.g. the Photonic Crystals book
by Joannopoulos et al.), so K perhaps isn't where you think it is.

In fact, K is the nearest-neighbor direction in the physical lattice
(since it is the longest vector in reciprocal space), which in the MPB
example is the +y and -y direction (along with all 60-degree rotations
thereof).

If you use the reciprocal->cartesian function to convert (-1/3, 1/3) to
Cartesian coordinates, you'll find that it corresponds to (0, -2/3), which
is one of the nearest-neighbor directions as expected.

Of course, there are also 5 other equivalent choices for K, and 12 choices
overall for the irreducible Brillouin zone.  This particular choice was
handed down to me by God (i.e. my advisor) via an input file for an older
band-structure code we were using when I first joined the MIT group in '96.

Cordially,
Steven G. Johnson

Another related post from 08/14/05 with subject Re: Question abt the triangular lattice model in MPB manual:

> 1. the reciprocal lattice for triangular lattice is still triangular lattice, which is rotated 30 degree about the center. So the angle between its two basis should not be 45 degree.

You are still thinking of a Cartesian coordinate system.  The angle between (0,0.5,0) and (-1/3,1/3,0) is *not* 45 degrees.  The reason is that the reciprocal-lattice-vector basis is not orthogonal.

For example, (-1/3,1/3,0) is -1/3 * G1 + 1/3 * G2 + 0 * G3 where G1,G2,G3 are the three reciprocal lattice vectors, which are *not* orthogonal axes.

> 2. the magnitude looks suspicious too: The coordinate for the first base based on the unit vector in k-space is (0 1/2 0) but the second one couldn't be (1/3 1/3 0) otherwise their magnitude are different, which should be identical for triangular lattice. The result I got was posted in my previous mail, with the title of correction.

Because it's not an orthogonal basis, the magnitude is not the square root of the sum of the squares of the components.  In any case, the M and K points are not *supposed* to have the same magnitude -- they are inequivalent points in the Brillouin zone.

Cordially,
Steven G. Johnson 
soamaven commented 5 years ago

Thanks @oskooi. The Nanocomp team is always so patient and helpful.

@JohnWeiner This link provides some further graphical interpretation of the (1/3, 1/3/ 0) K point.

this particular choice was handed down to me by God (i.e. my advisor)

:joy: This is how I feel being able to interact with Prof. Johnson through github, maillist. The God of CEM.

In fact, K is the nearest-neighbor direction in the physical lattice (since it is the longest vector in reciprocal space)

I think this answers my questions in the affirmative.

Because it's not an orthogonal basis, the magnitude is not the square root of the sum of the squares of the components.

This is something I once knew but am glad to be reminded of.