Eomys / pyleecan

Electrical engineering open-source software providing a user-friendly, unified, flexible simulation framework for the multiphysic design and optimization of electrical machines and drives
https://www.pyleecan.org
Apache License 2.0
151 stars 127 forks source link

[NF] Coupling with Elmer (Electromagnetics) #202

Open ajpina opened 3 years ago

ajpina commented 3 years ago

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

SebGue commented 3 years ago

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

BonneelP commented 3 years ago

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

SebGue commented 3 years ago

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

ajpina commented 3 years ago

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

SebGue commented 3 years ago

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

ajpina commented 3 years ago

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

BonneelP commented 3 years ago

Hello,

I took a look on your PR and it looks great :) Here are some modifications suggestions / remarks I've found to improve it:

Best regards, Pierre

ajpina commented 3 years ago

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

SebGue commented 3 years ago

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

BonneelP commented 3 years ago

Hello All,

+1 for the translation_dict :)

Anyway, as Sebastian suggested in #253 (and as we already talked in #195), we will set the standard for the labels once and for all and we can try to limit the labels to 30 characters.

Best regards, Pierre

SebGue commented 3 years ago

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

ajpina commented 3 years ago

Hi, I am working in mac os, I believe your (Windows) subprocess call doesn't find ElmerGrid in the PATH. I bet that if you add the whole absolute path to your ElmerGrid installation in the subprocess call (gen_elmer_mesh.py), it will work. I guess that probably your python console doesn't have the proper PATH where ElmerGrid can be found. Hope this helps. Best, Alejandro

On Thu, 3 Dec 2020, 04:22 Sebastian Günther, notifications@github.com wrote:

Hello,

@ajpina https://github.com/ajpina Could you tell me what your operating system is? @BonneelP https://github.com/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

— You are receiving this because you were mentioned. Reply to this email directly, view it on GitHub https://github.com/Eomys/pyleecan/issues/202#issuecomment-737775783, or unsubscribe https://github.com/notifications/unsubscribe-auth/AHYYMBOJ7W54PVNTBMBRBZDSS5KGVANCNFSM4TCJLNFQ .

SebGue commented 3 years ago

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

BonneelP commented 3 years ago

Hello Sebastian and Alejandro,

I run the test MagElmer on a Windows 10 - 64 bits and I got the WinError 2: image

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

SebGue commented 3 years ago

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

BonneelP commented 3 years ago

Hello Sebastian,

I'm on the master of pyleecan. Here is my Methods/Simulation/MagElmer/gen_elmer_mesh file before modification: image

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: image

Then when I run Tests/Functions/test_magelmer.py I got image

When I type "gmsh" in a Windows command the program opens normally. image

I can import the package (and last time I checked, the gmsh tests seems to run on my side) image

When I run your Test I got: image "File not found"

I can run other tests if you want.

Best regards, Pierre

SebGue commented 3 years ago

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

BonneelP commented 3 years ago

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

SebGue commented 3 years ago

Hello @ajpina

do you already have a parser for the 'SaveLine' results of Elmer?

Best regards, Sebastian

ajpina commented 3 years ago

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)

SebGue commented 3 years ago

Thank you, that's exactly what I was hoping for :-)

ajpina commented 3 years ago

Hi, Elmer developers changed their max number of characters for body names from 30 to 50 as you can read here: https://github.com/ElmerCSC/elmerfem/issues/249

Screen Shot 2020-12-08 at 4 25 16 PM (2) Screen Shot 2020-12-08 at 4 30 49 PM (2)

The picture shows the mesh in ElmerGUI after being converted from Gmsh format.

Best, Alejandro

SebGue commented 3 years ago

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

ajpina commented 3 years ago

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

SebGue commented 3 years ago

So we should not rely on this Elmer fix. I'm also not keen on compiling :-)

ajpina commented 3 years ago

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

EmileDvs commented 3 years ago

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