Making time-stepping FEM faster, part 4

Hi folks! After a couple of serious and not-so-serious lifestyle-esque posts, it’s time to get back to work. Mine, to be exact. If you’ve read the previous instalments in this series, you know that I’m currently working on a faster alternative to time-stepping FEM. If you haven’t, I suggest you start from the beginning. Today, it’s the time to discuss some implementation issues.

Previously

In the previous post, I described the basics of the harmonic balance method. The time-dependent vector potential $\avec(t)$ is written as a truncated Fourier series

$\avec(t) \approx \sum\limits_{k=1}^{2P+1} \avec_k \timefunction_k(t)$ ,

where $\timefunction_k (t)$ is a convenience notation for either $\cos$ or $\sin$ . This sum approximation is then substituted into the governing equation

$\stiffnessmatrix(t) \avec(t) = \loadvector(t)$ ,

and the unknown coefficients $\avec_k$ are then solved with the Galerkin’s method.

This gives us a big block matrix equation, with the matrix blocks $\stiffnessmatrix_{rc}$ equal to

$\stiffnessmatrix_{rc} = \langle \timefunction_r(t) \stiffnessmatrix(t) \timefunction_c(t) \rangle$ ,

where the angle brackets $\langle \rangle$ correspond to calculating the time-average. From this matrix system, the coefficients $\avec_k$ can then be solved.

Implementation

However, this leaves us with the not-so-trivial problem of actually computing these blocks. As I described earlier, the stiffness matrix $\stiffnessmatrix(t)$ can be divided into constant and time-dependent parts, i.e.

$\stiffnessmatrix(t) = \stiffnessmatrix_0 + \Sag(t)$ .

Obviously, the static component $\stiffnessmatrix_0$ is no problem to handle. However, this cannot be said for the time-dependent component $\Sag$ .

Indeed, the most typical approach to model the rotation of the machine is to change the finite element mesh in the air-gap, where the stator and rotor are moving past each other. In this approach (and most others, for that matter), the entries of $\Sag(t)$ will experience sharp, discontinuous changes as the elements of the mesh change. Thus, calculating the time-integral

$\langle \timefunction_r(t) \stiffnessmatrix(t) \timefunction_c(t) \rangle$

accurately could be quite complicated indeed.

However, and luckily, this is a problem I have already solved. Remember earlier when I spoke about calculating the Fourier series of the matrix $\Sag(t)$ ? In other words, I already know how to write $\Sag(t)$ as

$\Sag(t) \approx \sum\limits_{k=1}^{N} \Sag_k \timefunction_k(t)$ ,

and how to compute the coefficients $\Sag_k$ .

Number of required terms

Of course, having to compute a huge number $N$ of the terms would not be very practical, would it. I mean, that would be required to present the highly discontinuous entries of $\Sag$ accurately, right?

Well, yes and no. Accurate representation of $\Sag$ would indeed need quite a many terms. However, we are not interested in an accurate series approximation of $\Sag$ – we are interested in calculating the matrix blocks

$\langle \timefunction_r(t) \stiffnessmatrix(t) \timefunction_c(t) \rangle$ !

Thus, let’s now substitute the approximation for $\Sag(t)$ here, to get

$\langle \timefunction_r(t) \stiffnessmatrix(t) \timefunction_c(t) \rangle = \langle \timefunction_r(t) \timefunction_c(t) \rangle \stiffnessmatrix_0 + \sum\limits_{k=1}^N \langle \timefunction_r(t) \timefunction_k(t) \timefunction_c(t) \rangle \Sag_k$ .

In other words, we now longer have a time-average of a matrix – we have the time average of the product of three trigonometric functions $\langle \timefunction_r \timefunction_k \timefunction_c \rangle$ . And this can be expressed in closed form quite easily.

And even more better, there are only a limited number of terms $k$ for which the expression is non-zero. And these values of $k$ can also be determined analytically. Thus, it is usually sufficient to compute something like $N = 4P$ of the air-gap matrix Fourier coefficients.

And this is obviously a much easier task to computing a total of $(2P+1)^2$ of the matrix block integrals, like we would have been forced to do before.

Check out EMDtool - Electric Motor Design toolbox for Matlab.

Need help with electric motor design or design software? Let's get in touch - satisfaction guaranteed!

Making time-stepping FEM faster, part 4

3 thoughts on “Making time-stepping FEM faster, part 4”

Pavel Ponomarev says:

29th July 2016 at 11:32

Hi,

Great posts on Electrical Engineering! Those information bombs really help with understanding of those FEM tricks. Thank you! And please write more=)

BR, Pavel

1. Antti (admin) says:
  
  29th July 2016 at 13:02
  
  Hi Pavel,
  
  Thank you so much for your words! Nice that you like them 🙂
  
  And I will definitely continue writing 🙂
  
  -Antti
  
Antti (admin) says:

29th July 2016 at 13:03

And in general: if you have questions or ideas for blog topics, I’m always open for them 🙂

Antti Lehikoinen

Making time-stepping FEM faster, part 4

Previously

Implementation

Number of required terms

3 thoughts on “Making time-stepping FEM faster, part 4”

Leave a Reply Cancel reply