[NF] Coupling with Elmer (Electromagnetics)
ajpina opened this issue · 27 comments
Hi Team,
In order to develop a coupling with Elmer solver we need: GMSH -> ElmerGrid -> ElmerSolver -> VTK Postprocessing.
I can collaborate with the coupling for solving electromagnetics since GMSH coupling is already functional but It would be great if you create a skeleton for modules and functions so they follow the same principle as FEMM coupling. We can later work on the other Physics that can also be solved with Elmer.
Unfortunately, Elmer does not provide a python framework and the coupling needs to be performed through process calls. Please see an example on the following link: here
https://github.com/ajpina/emanfes/blob/master/emanfes/elmer/elmer_solver.py
I could populate the functions with similar code/approach once you define the placeholder for these according to pyleecan's module organization.
Best,
Alejandro
Hello Alejandro,
I already saw you got an Elmer solver in your other repo and wanted you to ask if you could contibute it to pyleecan as well. So that's great news.
I also had and discussion with Pierre recently about Elmer vs. GetDP. So why did you choose Elmer?
Last time i had a look at Elmer was maybe ten years ago :-) What are the current capabilities of Elmer EM right now? Is it AC or steady state, what element types and element orders are supported?
Do you also have other solves ready. I am particular interessted in thermal and structural simulations.
@BonneelP Do you think you could setup a skeleton beside your tutorials ;) ? Otherwise I can do it within the next days.
BTW: I think the Gmsh build dict and Elmer coupling should not be oriented too strictly to FEMM coupling since it will be a new class, e.g. MagElmer.
Best regards, Sebastian
Hello,
The coupling with Elmer would be a significant step for pyleecan :)
The webinars are taking longer than expected because we take the opportunity to improve and introduce some new feature alongside (spoiler, take a look at the GUI DXF_Hole ;) It is not finished yet). So I won't personally have the time to create the skeleton although I'm really interested to follow this issue. I should have the time to answer the question in this issue and maybe someone else at Eomys can take a look as well.
@SebGue, what do you have in mind for the skeleton ?
Best regards,
Pierre
Hello Pierre,
I think, I will first copy and paste MagFEMM, clear the files and add some comments (i.e. thoughts and ideas), so nothing special for now.
@ajpina What is the actual output of Elmer solver? Is there already some some post processing included, e.g. fluxlinkage calculation?
BR Sebastian
Hi,
Since the geometry is already built in Gmsh, you can take it forward to Getdp or Elmer for solving. I decided to go further with Elmer long time ago just because the documentation and the multi physic capability since you can invoke different solvers in cascade and they are fed by the results of the previous solutions.
It's my understanding that in terms of electromagnetics we can solve in Elmer so far, magneto static, transient and harmonic problems.
Most of my analysis have been with PM motors using transient solvers. You can also couple the electric circuit to it.
Since the geometry in draw_GMSH() can also output only lams, structural analysis can be performed with the same approach. Solvers for thermal problems are also available, so all these may come handy in the future.
I bet that Getdp can do similar jobs but I can't tell because I don't have the experience with it.
Regarding the outputs, there are several: Elmer provide some out of the box scalar outputs in a text file. Also, I used to include an arc in the air gap to sample the flux density and third, the results are written in VTU format, so you can access to magnetic vector potential or flux density at any time step.
I use VTK python libraries to extract MVP and calculate flux linkage from it. The VTU files are intended to be visualized in Paraview.
Hope this helps clear things out.
Thanks,
Alejandro
Hello Alejandro,
I created a PR with a skeleton for MagElmer class. I hope that it will fulfill your expectations somehow. If you need some further details or changes just let me know.
Best regards, Sebastian
Hi @SebGue , @BonneelP,
I just sent PR #236 with the first conversion from gmsh to elmer mesh via ElmerGrid. Could you guys create a method call gen_elmer_mesh(self, ouput) to the MagElmer object?
I would rather deal with Elmer mesh generation and Elmer settings in that method than including all that code in comp_flux_airgap().
Thanks,
Alejandro
Hello,
I took a look on your PR and it looks great :) Here are some modifications suggestions / remarks I've found to improve it:
- Tests/Functions/test_magelmer.py: when calling makedirs, always first use isdir, if the folder already exist it will create an error
- boundary_prop is very handy with the new label
- pyleecan/Functions/GMSH/draw_GMSH.py: We have a method "split_half" for all the Arc daughters that update the arc to either the first half (begin to middle) or the second half (middle to end). Maybe we can extend this method by adding a new parameter (return_both=False) to return both "half arc" instead of modifying the arc. (or maybe we can create a new method for that).
- pyleecan/Functions/GMSH/get_air_box.py line 28: What do you think of lam_ext = lam_list[-1] instead of lam_list[1] ? Do you think it would work with more than 2 laminations ?
- pyleecan/Methods/Simulation/MagElmer/solve_FEA.py: We have changed the way we handle the output in MagFEMM, now we define a dictionary with all the outputs, that way you can add Bz directly into the dictionary and create the data object in the store method of OutMag (sorry for changing again the way the models are organized. We are still experimenting to find the most convenient way to define several models that capitalize some code).
Best regards,
Pierre
Hi,
I am having issues when converting Gmsh to Elmer mesh because the surface labels are too long. ElmerGrid is able to translate names that are up to 30 characters and we have names that are way longer, particularly for the magnets and laminations, for instance: HoleMagnet_Rotor_Parallel_N_R0_T0_S0.
Now the problem is that I used to control the creation of ElmerSolver settings files by using the assigned names during Gmsh file generation, and now I have no way to control the random labels that Elmer assigns when the surfaces name's are longer than 30 characters and it creates a mismatch with those returned by build_geometry().
What do you guys suggest is the best course of actions?
Thanks,
Alejandro
Hello Alejandro,
what do you think of rename the surface labels in draw_GMSH just before setting the physical groups?
New names could simply be 'Magnet_1', 'Magnet_2', etc. Futher you could introduce a 'translation' dict to retrive the original names, e.g.
translation_dict = {'Magnet_1': 'HoleMagnet_Rotor_Parallel_N_R0_T0_S0'}
BTW: Good to know this limitation since I also started to utilize Elmer (and ElmerGrid).
Best Regards,
Sebastian
Hello,
@ajpina Could you tell me what your operating system is?
@BonneelP Could you run the MagElmer test and see if you also get an 'WinError 2 - file not found' error.
I am struggling with subprocess.Popen but don't know why. With some test script ElmerGrid is found but the actual implementation (in MagElmer and StructElmer) doesn't work although both seem to be absolutely the same.
Nevertheless, I think I already found a workaround...
Best regards, Sebastian
Hello,
I have Elmer in PATH, tried to move Elmer, set folder access permissions, ... nothing helped.
The strange thing is, that my test script is working while the actual implemetation doesn't. Although the code is exactly the same.
Also python can find ElmerGrid (i.e. distutils.spawn.find_executable method) and my windows shell find it too.
As a workaround, I did exactly what you suggest and now it works as expected.
Nevertheless its very strange and I am really interested if someone can reproduce this issue.
BR Sebastian
Hello Sebastian and Alejandro,
I run the test MagElmer on a Windows 10 - 64 bits and I got the WinError 2:
I have nothing (yet) related to Elmer installed on my computer but I have gmsh and the corresponding package installed.
Hope it helps you understand what is going on.
Best regards,
Pierre
Hello Pierre,
unfortunatly no. Maybe you could do the following, which should be equivalent:
Temporarly replace the 'cmd_elmergrid' in gen_elmer_mesh method (around line 44) and run the test.
from distutils import spawn
cmd_elmergrid = ["gmsh"]
binary = spawn.find_executable(cmd_elmergrid[0])
if not binary:
raise Exception
This should raise an exception if GMSH is not in PATH (this should not happen), raise the WinError2 (this would reproduce my error) or raise an other error due to missing mesh files (i.e. POpen is executed as expected).
Further you could run the following code:
import subprocess
from os import getcwd
from distutils import spawn
class TestElmerGrid:
def create_mesh(self):
cmd = ["gmsh"]
process = subprocess.Popen(
cmd,
stdout=subprocess.PIPE,
stderr=subprocess.PIPE, # shell=True
)
(stdout, stderr) = process.communicate()
process.wait()
if __name__ == "__main__":
obj = TestElmerGrid()
obj.create_mesh()
This code runs as expected on my system and that's the strange thing, since it is the same code snippet as in gen_elmer_mesh.
My actual PR fix this issue, but I really would like to know why the above issue is going to happen.
Best regards, Sebastian
Hello Sebastian,
I'm on the master of pyleecan. Here is my Methods/Simulation/MagElmer/gen_elmer_mesh file before modification:
cmd_elmergrid is on line 30 for me. So if I understood well (which maybe not the case) here is how I modified the script:
Then when I run Tests/Functions/test_magelmer.py I got
When I type "gmsh" in a Windows command the program opens normally.
I can import the package (and last time I checked, the gmsh tests seems to run on my side)
When I run your Test I got:
"File not found"
I can run other tests if you want.
Best regards,
Pierre
Oh, I forgot... Did you import distutils spawn? But I guess you did.
So it seems your python doesn't even look at your PATH (because of fig 3 exeception).
Maybe we shouldn't investigate further not to waste too much time since I have a fix for this issue?!
Best regards, Sebastian
Yes the "from distutils import spawn" was in my script but higher than my screenshot ;)
If you already have found a solution, I think we should focus on providing unittest for this solution to make sure that it works.
We were thinking lately about adding a method to check that everything is correctly set for each model (something like "check_input" to call at the beginning of run) maybe we can investigate how we can check that everything is installed as it should (gmsh/elmer...).
This check_input method would also take all the verifications that we currently do in the gen_input and adapt it to the current model (since not all models have the same requirements).
Best regards,
Pierre
Hello @ajpina
do you already have a parser for the 'SaveLine' results of Elmer?
Best regards, Sebastian
I have this, it should do the work:
ecp, mfe, agt, iv, im, tq = np.loadtxt(scalars_file, unpack=True, usecols=(0, 1, 2, 3, 4, 5))
# This order must match scalars.dat.names
# 1: res: eddy current power
# 2: res: magnetic field energy
# 3: res: air gap torque
# 4: res: inertial volume
# 5: res: inertial moment
# 6: res: group 1 torque
cls.time_axis = np.linspace(0, time_step * (steps - skip_steps - 1), steps - skip_steps)
cls.angle_axis = cls.time_axis * degrees_per_sec
cls.torque_elmer = agt[skip_steps:] * stack_length
cls.torque_akkio = tq[skip_steps:] * fractions * stack_length
data = np.loadtxt(line_file, usecols=(0, 4, 5, 7, 8, 10, 11))
# This order must match lines.dat.names
# 1: Time step
# 2: Iteration step
# 3: Boundary condition
# 4: Node index
# 5: coordinate 1
# 6: coordinate 2
# 7: coordinate 3
# 8: magnetic flux density 1
# 9: magnetic flux density 2
# 10: magnetic flux density 3
# 11: magnetic flux density e 1
# 12: magnetic flux density e 2
# 13: magnetic flux density e 3
Br_list = []
Bt_list = []
theta_list = []
for i in range(skip_steps, steps + 1):
index = np.where(data == i)
data_step = data[index[0]]
rho = np.sqrt(data_step[:, 1] ** 2 + data_step[:, 2] ** 2)
theta = np.arctan2(data_step[:, 2], data_step[:, 1])
Br1 = np.multiply(data_step[:, 3], np.cos(theta)) + np.multiply(data_step[:, 4], np.sin(theta))
Bt1 = -np.multiply(data_step[:, 3], np.sin(theta)) + np.multiply(data_step[:, 4], np.cos(theta))
max_theta = np.max(theta)
min_theta = np.min(theta)
data_step_polar = np.stack((rho, theta, Br1, Bt1), axis=-1)
data_step_polar = data_step_polar[data_step_polar[:, 1].argsort()]
cs_r = CubicSpline(data_step_polar[:, 1], data_step_polar[:, 2])
cs_t = CubicSpline(data_step_polar[:, 1], data_step_polar[:, 3])
theta_fine = np.linspace(min_theta, max_theta, num=int(2880 / fractions))
Br_list.append(cs_r(theta_fine))
Bt_list.append(cs_t(theta_fine))
theta_list.append(theta_fine)
Thank you, that's exactly what I was hoping for :-)
Hi,
Elmer developers changed their max number of characters for body names from 30 to 50 as you can read here:
ElmerCSC/elmerfem#249
The picture shows the mesh in ElmerGUI after being converted from Gmsh format.
Best,
Alejandro
Hi Alejandro,
I already saw your issue at Elmer repo. Thanks for noticing.
Do you have any idea when they will update their binaries? Do you compile yours yourself?
Best, Sebastian
It will take a while. Elmer developers take several months, sometime years, to release new binaries (versions). I am building from the source (potential v9) which is not fun at all in mac os (too many problems with dependencies). Pretty sure that compiling it in Linux or Windows would be much easier.
Best,
Alejandro
So we should not rely on this Elmer fix. I'm also not keen on compiling :-)
Hi @BonneelP , @EmileDvs ,
How can I get the stator flux A-axis? I found these methods stator.comp_angle_d_axis() , stator.comp_angle_q_axis() and machine.comp_angle_offset_initial() but their results don't make sense to me.
Could you guys advice how to calculate the rotor initial position so that the rotor d-axis is aligned with stator A-axis?
Thanks,
Alejandro
Hi @ajpina,
The methods you found are the ones you need to use indeed. machine.comp_angle_offset_initial() returns the angle shift between rotor d-axis and stator d-axis (same as A-axis for stator as it doesn't rotate) so that they are well aligned in MagFEMM simulation for example.
Please let me know if you need further information on this,
Emile