nanoHUB U Nanophotonic Modeling/Nanophotonic Modeling Lecture 4.9: Improving FEM Accuracy
========================================
 

[Slide 1] Hey, everyone. Welcome to Lecture 4.9. In the last few lectures, we kind of developed the finite element method and showed how we could use the Galerkin Method to achieve high accuracies in a systematic fashion. And so now we're actually going to explore how does this connect with optics. 

[Slide 2] So of course one of the key challenges in optics is that any kind of fields that we're generating can potentially be reflected back into our spatial grid. And that's true whether it's finite difference time domain or if it's finite element method time domain problem. And so in order to prevent their reflections then we actually need to add a new feature to the finite element method which we can call perfectly matched layers. And so much like we saw in finite difference time domain we can introduce this by basically a few different methods. But we're going to use a slightly different approach than we had for finite difference which is to actually stretch the coordinates of the system. And so this basically is effected by transforming the gradient operator into A dotted with the gradient operator. So we can represent our A coefficient as basically having a small imaginary component. So that means we have both real and imaginary space at this point. And this beta actually is chosen very specifically based on the details of our simulation, and so this is kind of the formula that we typically use. This is explained in detail in Obaya's book, if you want to look it up. 

[Slide 3] And so then in the context of finite element propagation method, then of course, remember that we already had written down a solution. And we had also written down the expressions for propagation, and we can even now include a second derivative of phi which looks like this okay? So then this actually is going to be a little more accurate than the earlier versions of finite method that we used which will only have the first order derivative. And in order to get a reasonable solution, which accounts for both perfectly matched layers and finite element methods with second derivatives, then we'll use what's called a Pade approximant. Which would be based on this initial value of this derivative of the envelope function, which we call zeta. Okay, so this value here, zeta, which would be dphi/dz. So, if we choose that the initial second derivative of the envelope function would be a fixed value then presumably zeta of 0 would be 0. And so then we can actually analytically calculate that zeta should look something like this, which J beta times square root of these terms divided by beta squared minus one times the value of phi. And this is kind of related to applying a formal solution of phi to zeta. 

[Slide 4] Okay, and then we can go back to what we had talked about in the last couple of lectures to apply the Galerkin method. And then we can actually write down a set of values for the fields. Here we will assume that the fields are H fields. And then we have basically a set of coefficients in the z direction times some transverse coordinates phi that look like 5j of xy times some unit vectors in the x and y directions. And then that basically allows us to rewrite the overall set of equations in beam propagation method. As mass matrix times the second derivative of hT in the z direction, minus 2 gamma, times another mass matrix. Again, times the first derivative of hT in the z direction. And then plus kind of a spring constant type equation, plus gamma squared, times M. So this is actually very analogous in many ways to what we saw in the case of the harmonic oscillator, except we have this extra mass matrix, and then we set that equal to zero. Of course, it won't necessarily be exactly zero, but then we can define that as error function, and then we can integrate with respect to the test cases. Here It's showing the details of what the mass matrix looks like. So again, it's like the overlap of the different basis elements. And then k, i, j, is basically looking at what happens when you take the curl of those elements and then subtract them. So, obviously, the mathematics is a bit complex but you don't have to memorize the details. 

[Slide 5] But the basic point is that we have a fairly reasonable method at this point to calculate. Find finite element method, beam propagation method. And we probably are most concerned about, how accurate is this approach? So, of course it may depend a lot on the basis functions. But if we assume that we have good basis functions, then we have a few options to improve the performance. One of course would be just decrease the mesh size, which we showed in the last lecture. It's oftentimes a very viable option, and can systematically reduce the errors basically to zero. That's called H-adaptivity. There's also the so called P-adaptivity which allows you to increase the degree of the fitted polynomials, often times starting with linear basis functions. But if there's high curvature in certain regions then you might want to go to quadratic or even quartic terms. But then of course you can combine both of these types, of adaptivity into what's called HP-adaptivity, in order to get even lower errors. 

[Slide 6] And, so, typically the strategy would be to create an initial mesh, and then to compute your first solution on that mesh which may be very coarse and then compute the error and of course accept that the error may be much larger than you would desire as a final result. So, then you can basically check how well it matches up and then refine the mesh spacing if needed using the H-abductivity method or if you're using HP you can simultaneously both increase the polynomial degree and the mesh. Now of course the catch with the HP method is it can be very computationally expensive because the calculation of each of the polynomial elements becomes more expensive as you add more coefficients. And at the same time you've added more elements. So the time could increase very quickly. In the next lecture, we'll talk about more details of the finite element method.