From 3b653b78f96dfaed61d78db31b52a226983f0782 Mon Sep 17 00:00:00 2001 From: Uwe Date: Sun, 5 Feb 2023 07:46:12 +0100 Subject: [PATCH] [FEM] Elmer: add support for 3D magnetodynamics - adds the corresponding Elmer equation (it is now possible to do Elmer's tutorial example no. 14) --- src/Mod/Fem/CMakeLists.txt | 2 + src/Mod/Fem/Gui/Command.cpp | 26 +- src/Mod/Fem/Gui/Resources/Fem.qrc | 1 + .../icons/FEM_EquationMagnetodynamic.svg | 97 +++++ .../Fem/Gui/Resources/ui/CurrentDensity.ui | 10 +- .../Resources/ui/ElectrostaticPotential.ui | 18 +- src/Mod/Fem/ObjectsFem.py | 14 + src/Mod/Fem/femcommands/commands.py | 21 + .../elmer/equations/magnetodynamic.py | 206 ++++++++++ .../elmer/equations/magnetodynamic_writer.py | 375 ++++++++++++++++++ src/Mod/Fem/femsolver/elmer/writer.py | 26 ++ src/Mod/Fem/femsolver/equationbase.py | 10 + src/Mod/Fem/femtest/app/test_object.py | 38 +- src/Mod/Fem/femtest/app/test_open.py | 5 + src/Mod/Fem/femtest/gui/test_open.py | 5 + 15 files changed, 837 insertions(+), 17 deletions(-) create mode 100644 src/Mod/Fem/Gui/Resources/icons/FEM_EquationMagnetodynamic.svg create mode 100644 src/Mod/Fem/femsolver/elmer/equations/magnetodynamic.py create mode 100644 src/Mod/Fem/femsolver/elmer/equations/magnetodynamic_writer.py diff --git a/src/Mod/Fem/CMakeLists.txt b/src/Mod/Fem/CMakeLists.txt index 65ce6b8095..6ed260d2df 100755 --- a/src/Mod/Fem/CMakeLists.txt +++ b/src/Mod/Fem/CMakeLists.txt @@ -266,6 +266,8 @@ SET(FemSolverElmerEquations_SRCS femsolver/elmer/equations/heat.py femsolver/elmer/equations/heat_writer.py femsolver/elmer/equations/linear.py + femsolver/elmer/equations/magnetodynamic.py + femsolver/elmer/equations/magnetodynamic_writer.py femsolver/elmer/equations/magnetodynamic2D.py femsolver/elmer/equations/magnetodynamic2D_writer.py femsolver/elmer/equations/nonlinear.py diff --git a/src/Mod/Fem/Gui/Command.cpp b/src/Mod/Fem/Gui/Command.cpp index 2ff8be389e..15e632524e 100644 --- a/src/Mod/Fem/Gui/Command.cpp +++ b/src/Mod/Fem/Gui/Command.cpp @@ -1386,6 +1386,8 @@ void CmdFEMCompEmEquations::activated(int iMsg) else if (iMsg == 1) rcCmdMgr.runCommandByName("FEM_EquationElectricforce"); else if (iMsg == 2) + rcCmdMgr.runCommandByName("FEM_EquationMagnetodynamic"); + else if (iMsg == 3) rcCmdMgr.runCommandByName("FEM_EquationMagnetodynamic2D"); else return; @@ -1410,7 +1412,9 @@ Gui::Action* CmdFEMCompEmEquations::createAction() QAction* cmd1 = pcAction->addAction(QString()); cmd1->setIcon(Gui::BitmapFactory().iconFromTheme("FEM_EquationElectricforce")); QAction* cmd2 = pcAction->addAction(QString()); - cmd2->setIcon(Gui::BitmapFactory().iconFromTheme("FEM_EquationMagnetodynamic2D")); + cmd2->setIcon(Gui::BitmapFactory().iconFromTheme("FEM_EquationMagnetodynamic")); + QAction* cmd3 = pcAction->addAction(QString()); + cmd3->setIcon(Gui::BitmapFactory().iconFromTheme("FEM_EquationMagnetodynamic2D")); _pcAction = pcAction; languageChange(); @@ -1456,15 +1460,27 @@ void CmdFEMCompEmEquations::languageChange() EquationElectricforce->getStatusTip())); } + Gui::Command* EquationMagnetodynamic = + rcCmdMgr.getCommandByName("FEM_EquationMagnetodynamic"); + if (EquationMagnetodynamic) { + QAction* cmd2 = a[2]; + cmd2->setText(QApplication::translate("FEM_EquationMagnetodynamic", + EquationMagnetodynamic->getMenuText())); + cmd2->setToolTip(QApplication::translate("FEM_EquationMagnetodynamic", + EquationMagnetodynamic->getToolTipText())); + cmd2->setStatusTip(QApplication::translate("FEM_EquationMagnetodynamic", + EquationMagnetodynamic->getStatusTip())); + } + Gui::Command* EquationMagnetodynamic2D = rcCmdMgr.getCommandByName("FEM_EquationMagnetodynamic2D"); if (EquationMagnetodynamic2D) { - QAction* cmd2 = a[2]; - cmd2->setText(QApplication::translate("FEM_EquationMagnetodynamic2D", + QAction* cmd3 = a[3]; + cmd3->setText(QApplication::translate("FEM_EquationMagnetodynamic2D", EquationMagnetodynamic2D->getMenuText())); - cmd2->setToolTip(QApplication::translate("FEM_EquationMagnetodynamic2D", + cmd3->setToolTip(QApplication::translate("FEM_EquationMagnetodynamic2D", EquationMagnetodynamic2D->getToolTipText())); - cmd2->setStatusTip(QApplication::translate("FEM_EquationMagnetodynamic2D", + cmd3->setStatusTip(QApplication::translate("FEM_EquationMagnetodynamic2D", EquationMagnetodynamic2D->getStatusTip())); } } diff --git a/src/Mod/Fem/Gui/Resources/Fem.qrc b/src/Mod/Fem/Gui/Resources/Fem.qrc index 3636ce6723..7f7fdcc322 100755 --- a/src/Mod/Fem/Gui/Resources/Fem.qrc +++ b/src/Mod/Fem/Gui/Resources/Fem.qrc @@ -47,6 +47,7 @@ icons/FEM_EquationFlow.svg icons/FEM_EquationFlux.svg icons/FEM_EquationHeat.svg + icons/FEM_EquationMagnetodynamic.svg icons/FEM_EquationMagnetodynamic2D.svg diff --git a/src/Mod/Fem/Gui/Resources/icons/FEM_EquationMagnetodynamic.svg b/src/Mod/Fem/Gui/Resources/icons/FEM_EquationMagnetodynamic.svg new file mode 100644 index 0000000000..1e1b76239d --- /dev/null +++ b/src/Mod/Fem/Gui/Resources/icons/FEM_EquationMagnetodynamic.svg @@ -0,0 +1,97 @@ + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + image/svg+xml + + + + [Alexander Gryson] + + + 2017-03-11 + http://www.freecadweb.org/wiki/index.php?title=Artwork + + + FreeCAD + + + FreeCAD/src/Mod/ + + + FreeCAD LGPL2+ + + + https://www.gnu.org/copyleft/lesser.html + + + [agryson] Alexander Gryson + + + + + + + + + + + + + + + + + + + + + + + + + diff --git a/src/Mod/Fem/Gui/Resources/ui/CurrentDensity.ui b/src/Mod/Fem/Gui/Resources/ui/CurrentDensity.ui index 4613750d18..f5ea8fa32e 100644 --- a/src/Mod/Fem/Gui/Resources/ui/CurrentDensity.ui +++ b/src/Mod/Fem/Gui/Resources/ui/CurrentDensity.ui @@ -63,7 +63,10 @@ - Real part of potential x-component + Real part of potential x-component +Note: if a face was selected this will be the value + in normal face direction + settings for y and z will be ignored true @@ -107,7 +110,10 @@ - Imaginary part of potential x-component + Imaginary part of potential x-component +Note: if a face was selected this will be the value + in normal face direction + settings for y and z will be ignored true diff --git a/src/Mod/Fem/Gui/Resources/ui/ElectrostaticPotential.ui b/src/Mod/Fem/Gui/Resources/ui/ElectrostaticPotential.ui index ca783344d6..056bfeddcc 100644 --- a/src/Mod/Fem/Gui/Resources/ui/ElectrostaticPotential.ui +++ b/src/Mod/Fem/Gui/Resources/ui/ElectrostaticPotential.ui @@ -241,7 +241,8 @@ with a harmonic/oscillating driving force - Real part of potential x-component + Real part of potential x-component +Note: has no effect if a solid was selected true @@ -285,7 +286,8 @@ with a harmonic/oscillating driving force - Imaginary part of potential x-component + Imaginary part of potential x-component +Note: has no effect if a solid was selected true @@ -339,7 +341,8 @@ with a harmonic/oscillating driving force - Real part of potential y-component + Real part of potential y-component +Note: has no effect if a solid was selected true @@ -383,7 +386,8 @@ with a harmonic/oscillating driving force - Imaginary part of potential y-component + Imaginary part of potential y-component +Note: has no effect if a solid was selected true @@ -437,7 +441,8 @@ with a harmonic/oscillating driving force - Real part of potential z-component + Real part of potential z-component +Note: has no effect if a solid was selected true @@ -481,7 +486,8 @@ with a harmonic/oscillating driving force - Imaginary part of potential z-component + Imaginary part of potential z-component +Note: has no effect if a solid was selected true diff --git a/src/Mod/Fem/ObjectsFem.py b/src/Mod/Fem/ObjectsFem.py index e6d00b2c31..ec7d4f2f10 100644 --- a/src/Mod/Fem/ObjectsFem.py +++ b/src/Mod/Fem/ObjectsFem.py @@ -817,6 +817,20 @@ def makeEquationHeat( return obj +def makeEquationMagnetodynamic( + doc, + base_solver=None, + name="Magnetodynamic" +): + """makeEquationMagnetodynamic(document, [base_solver], [name]): + creates a FEM magnetodynamic equation for a solver""" + from femsolver.elmer.equations import magnetodynamic + obj = magnetodynamic.create(doc, name) + if base_solver: + base_solver.addObject(obj) + return obj + + def makeEquationMagnetodynamic2D( doc, base_solver=None, diff --git a/src/Mod/Fem/femcommands/commands.py b/src/Mod/Fem/femcommands/commands.py index 4e625983db..ed1f4e5d1b 100644 --- a/src/Mod/Fem/femcommands/commands.py +++ b/src/Mod/Fem/femcommands/commands.py @@ -533,6 +533,23 @@ class _EquationHeat(CommandManager): self.do_activated = "add_obj_on_gui_selobj_noset_edit" +class _EquationMagnetodynamic(CommandManager): + "The FEM_EquationMagnetodynamic command definition" + + def __init__(self): + super(_EquationMagnetodynamic, self).__init__() + self.menutext = Qt.QT_TRANSLATE_NOOP( + "FEM_EquationMagnetodynamic", + "Magnetodynamic equation" + ) + self.tooltip = Qt.QT_TRANSLATE_NOOP( + "FEM_EquationMagnetodynamic", + "Creates a FEM equation for\n magentodynamic forces" + ) + self.is_active = "with_solver_elmer" + self.do_activated = "add_obj_on_gui_selobj_noset_edit" + + class _EquationMagnetodynamic2D(CommandManager): "The FEM_EquationMagnetodynamic2D command definition" @@ -1239,6 +1256,10 @@ FreeCADGui.addCommand( "FEM_EquationHeat", _EquationHeat() ) +FreeCADGui.addCommand( + "FEM_EquationMagnetodynamic", + _EquationMagnetodynamic() +) FreeCADGui.addCommand( "FEM_EquationMagnetodynamic2D", _EquationMagnetodynamic2D() diff --git a/src/Mod/Fem/femsolver/elmer/equations/magnetodynamic.py b/src/Mod/Fem/femsolver/elmer/equations/magnetodynamic.py new file mode 100644 index 0000000000..e5a1aed498 --- /dev/null +++ b/src/Mod/Fem/femsolver/elmer/equations/magnetodynamic.py @@ -0,0 +1,206 @@ +# *************************************************************************** +# * Copyright (c) 2023 Uwe Stöhr * +# * * +# * This file is part of the FreeCAD CAx development system. * +# * * +# * This program is free software; you can redistribute it and/or modify * +# * it under the terms of the GNU Lesser General Public License (LGPL) * +# * as published by the Free Software Foundation; either version 2 of * +# * the License, or (at your option) any later version. * +# * for detail see the LICENCE text file. * +# * * +# * This program is distributed in the hope that it will be useful, * +# * but WITHOUT ANY WARRANTY; without even the implied warranty of * +# * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the * +# * GNU Library General Public License for more details. * +# * * +# * You should have received a copy of the GNU Library General Public * +# * License along with this program; if not, write to the Free Software * +# * Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 * +# * USA * +# * * +# *************************************************************************** + +__title__ = "FreeCAD FEM solver Elmer equation object Magnetodynamic" +__author__ = "Uwe Stöhr" +__url__ = "https://www.freecadweb.org" +# +## \addtogroup FEM +# @{ + +from femtools import femutils +from . import nonlinear +from ... import equationbase + +def create(doc, name="Magnetodynamic"): + return femutils.createObject( + doc, name, Proxy, ViewProxy) + + +class Proxy(nonlinear.Proxy, equationbase.MagnetodynamicProxy): + + Type = "Fem::EquationElmerMagnetodynamic" + + def __init__(self, obj): + super(Proxy, self).__init__(obj) + + obj.addProperty( + "App::PropertyBool", + "IsHarmonic", + "Magnetodynamic", + "If the magnetic source is harmonically driven" + ) + obj.addProperty( + "App::PropertyFrequency", + "AngularFrequency", + "Magnetodynamic", + "Frequency of the driving current" + ) + obj.addProperty( + "App::PropertyBool", + "UsePiolaTransform", + "Magnetodynamic", + "Must be True if basis functions for edge element interpolation\n" + "are selected to be members of optimal edge element family\n" + "or if second-order approximation is used." + ) + obj.addProperty( + "App::PropertyBool", + "QuadraticApproximation", + "Magnetodynamic", + "Enables second-order approximation of driving current" + ) + obj.addProperty( + "App::PropertyBool", + "StaticConductivity", + "Magnetodynamic", + "See Elmer models manual for info" + ) + obj.addProperty( + "App::PropertyBool", + "FixInputCurrentDensity", + "Magnetodynamic", + "Ensures divergence-freeness of current density" + ) + obj.addProperty( + "App::PropertyBool", + "AutomatedSourceProjectionBCs", + "Magnetodynamic", + "See Elmer models manual for info" + ) + obj.addProperty( + "App::PropertyBool", + "UseLagrangeGauge", + "Magnetodynamic", + "See Elmer models manual for info" + ) + obj.addProperty( + "App::PropertyFloat", + "LagrangeGaugePenalizationCoefficient", + "Magnetodynamic", + "See Elmer models manual for info" + ) + obj.addProperty( + "App::PropertyBool", + "UseTreeGauge", + "Magnetodynamic", + "See Elmer models manual for info\n" + "Will be ignored if 'UsePiolaTransform' is True" + ) + obj.addProperty( + "App::PropertyBool", + "LinearSystemRefactorize", + "Linear System", + "" + ) + + obj.IsHarmonic = False + obj.AngularFrequency = 0 + obj.Priority = 10 + + # the post processor options + obj.addProperty( + "App::PropertyBool", + "CalculateCurrentDensity", + "Results", + "" + ) + obj.addProperty( + "App::PropertyBool", + "CalculateElectricField", + "Results", + "" + ) + obj.addProperty( + "App::PropertyBool", + "CalculateElementalFields", + "Results", + "" + ) + obj.addProperty( + "App::PropertyBool", + "CalculateHarmonicLoss", + "Results", + "" + ) + obj.addProperty( + "App::PropertyBool", + "CalculateJouleHeating", + "Results", + "" + ) + obj.addProperty( + "App::PropertyBool", + "CalculateMagneticFieldStrength", + "Results", + "" + ) + obj.addProperty( + "App::PropertyBool", + "CalculateMaxwellStress", + "Results", + "" + ) + obj.addProperty( + "App::PropertyBool", + "CalculateNodalFields", + "Results", + "" + ) + obj.addProperty( + "App::PropertyBool", + "CalculateNodalForces", + "Results", + "" + ) + obj.addProperty( + "App::PropertyBool", + "CalculateNodalHeating", + "Results", + "" + ) + obj.addProperty( + "App::PropertyBool", + "DiscontinuousBodies", + "Results", + "" + ) + obj.CalculateCurrentDensity = False + obj.CalculateElectricField = False + # FIXME: at the moment FreeCAD's post processor cannot display elementary field + # results, therefore disable despite this is by default on in Elmer + obj.CalculateElementalFields = False + obj.CalculateHarmonicLoss = False + obj.CalculateJouleHeating = False + obj.CalculateMagneticFieldStrength = False + obj.CalculateMaxwellStress = False + obj.CalculateNodalFields = True + obj.CalculateNodalForces = False + obj.CalculateNodalHeating = False + obj.DiscontinuousBodies = False + + +class ViewProxy(nonlinear.ViewProxy, equationbase.MagnetodynamicViewProxy): + pass + +## @} diff --git a/src/Mod/Fem/femsolver/elmer/equations/magnetodynamic_writer.py b/src/Mod/Fem/femsolver/elmer/equations/magnetodynamic_writer.py new file mode 100644 index 0000000000..c8e74fdc92 --- /dev/null +++ b/src/Mod/Fem/femsolver/elmer/equations/magnetodynamic_writer.py @@ -0,0 +1,375 @@ +# *************************************************************************** +# * Copyright (c) 2023 Uwe Stöhr * +# * * +# * This file is part of the FreeCAD CAx development system. * +# * * +# * This program is free software; you can redistribute it and/or modify * +# * it under the terms of the GNU Lesser General Public License (LGPL) * +# * as published by the Free Software Foundation; either version 2 of * +# * the License, or (at your option) any later version. * +# * for detail see the LICENCE text file. * +# * * +# * This program is distributed in the hope that it will be useful, * +# * but WITHOUT ANY WARRANTY; without even the implied warranty of * +# * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the * +# * GNU Library General Public License for more details. * +# * * +# * You should have received a copy of the GNU Library General Public * +# * License along with this program; if not, write to the Free Software * +# * Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 * +# * USA * +# * * +# *************************************************************************** + +__title__ = "FreeCAD FEM Magnetodynamics Elmer writer" +__author__ = "Uwe Stöhr" +__url__ = "https://www.freecad.org" + +## \addtogroup FEM +# @{ + +from FreeCAD import Console +from FreeCAD import Units + +from .. import sifio +from .. import writer as general_writer + +class MgDynwriter: + + def __init__(self, writer, solver): + self.write = writer + self.solver = solver + + def getMagnetodynamicSolver(self, equation): + # output the equation parameters + s = self.write.createNonlinearSolver(equation) + if not equation.IsHarmonic: + s["Equation"] = "MgDyn" + s["Procedure"] = sifio.FileAttr("MagnetoDynamics/WhitneyAVSolver") + s["Variable"] = "av" + else: + s["Equation"] = "MgDynHarmonic" + s["Procedure"] = sifio.FileAttr("MagnetoDynamics/WhitneyAVHarmonicSolver") + s["Variable"] = "av[av re:1 av im:1]" + # round to get rid of numerical artifacts + frequency = float(Units.Quantity(equation.AngularFrequency).Value) + s["Angular Frequency"] = round(frequency, 6) + s["Exec Solver"] = "Always" + s["Optimize Bandwidth"] = True + s["Stabilize"] = equation.Stabilize + if equation.LinearSystemRefactorize is True: + s["Linear System Refactorize"] = True + if equation.UsePiolaTransform is True: + s["Use Piola Transform"] = True + if equation.QuadraticApproximation is True: + s["Quadratic Approximation"] = True + if equation.StaticConductivity is True: + s["Static Conductivity"] = True + if equation.FixInputCurrentDensity is True: + s["Fix Input Current Density"] = True + if equation.AutomatedSourceProjectionBCs is True: + s["Automated Source Projection BCs"] = True + if equation.UseLagrangeGauge is True: + s["Use Lagrange Gauge"] = True + if equation.LagrangeGaugePenalizationCoefficient != 0.0: + s["Lagrange Gauge Penalization Coefficient"] = \ + equation.LagrangeGaugePenalizationCoefficient + if equation.UseTreeGauge is True: + s["Use Tree Gauge"] = True + return s + + def getMagnetodynamicSolverPost(self, equation): + # output the equation parameters + s = self.write.createNonlinearSolver(equation) + s["Equation"] = "MgDynPost" + s["Exec Solver"] = "Before Saving" + s["Procedure"] = sifio.FileAttr("MagnetoDynamics/MagnetoDynamicsCalcFields") + if equation.IsHarmonic: + frequency = float(Units.Quantity(equation.AngularFrequency).Value) + s["Angular Frequency"] = round(frequency, 6) + s["Potential Variable"] = "av" + if equation.CalculateCurrentDensity is True: + s["Calculate Current Density"] = True + if equation.CalculateElectricField is True: + s["Calculate Electric Field"] = True + if equation.CalculateElementalFields is False: + s["Calculate Elemental Fields"] = False + if equation.CalculateHarmonicLoss is True: + s["Calculate Harmonic Loss"] = True + if equation.CalculateJouleHeating is True: + s["Calculate Joule Heating"] = True + if equation.CalculateMagneticFieldStrength is True: + s["Calculate Magnetic Field Strength"] = True + if equation.CalculateMaxwellStress is True: + s["Calculate Maxwell Stress"] = True + if equation.CalculateNodalFields is False: + s["Calculate Nodal Fields"] = False + if equation.CalculateNodalForces is True: + s["Calculate Nodal Forces"] = True + if equation.CalculateNodalHeating is True: + s["Calculate Nodal Heating"] = True + if equation.DiscontinuousBodies is True: + s["Discontinuous Bodies"] = True + s["Optimize Bandwidth"] = True + s["Stabilize"] = equation.Stabilize + return s + + def handleMagnetodynamicConstants(self): + permeability = self.write.convert( + self.write.constsdef["PermeabilityOfVacuum"], + "M*L/(T^2*I^2)" + ) + # we round in the following to get rid of numerical artifacts + self.write.constant("Permeability Of Vacuum", round(permeability, 20)) + + permittivity = self.write.convert( + self.write.constsdef["PermittivityOfVacuum"], + "T^4*I^2/(L^3*M)" + ) + self.write.constant("Permittivity Of Vacuum", round(permittivity, 20)) + + def handleMagnetodynamicMaterial(self, bodies): + # check that all bodies have a set material + for name in bodies: + if self.write.getBodyMaterial(name) == None: + raise general_writer.WriteError( + "The body {} is not referenced in any material.\n\n".format(name) + ) + for obj in self.write.getMember("App::MaterialObject"): + m = obj.Material + refs = ( + obj.References[0][1] + if obj.References + else self.write.getAllBodies()) + for name in (n for n in refs if n in bodies): + if "ElectricalConductivity" not in m: + Console.PrintMessage("m: {}\n".format(m)) + raise general_writer.WriteError( + "The electrical conductivity must be specified for all materials.\n\n" + ) + if "RelativePermeability" not in m: + Console.PrintMessage("m: {}\n".format(m)) + raise general_writer.WriteError( + "The relative permeability must be specified for all materials.\n\n" + ) + self.write.material(name, "Name", m["Name"]) + conductivity = self.write.convert(m["ElectricalConductivity"], "T^3*I^2/(L^3*M)") + conductivity = round(conductivity, 10) # to get rid of numerical artifacts + self.write.material(name, "Electric Conductivity", conductivity) + self.write.material( + name, "Relative Permeability", + float(m["RelativePermeability"]) + ) + # permittivity might be necessary for the post processor + if "RelativePermittivity" in m: + self.write.material( + name, "Relative Permittivity", + float(m["RelativePermittivity"]) + ) + + def _outputMagnetodynamicBodyForce(self, obj, name, equation): + if hasattr(obj, "CurrentDensity_re_1"): + # output only if current density is enabled and needed + if not obj.CurrentDensity_re_1_Disabled: + currentDensity = float(obj.CurrentDensity_re_1.getValueAs("A/m^2")) + self.write.bodyForce(name, "Current Density 1", round(currentDensity, 6)) + if not obj.CurrentDensity_re_2_Disabled: + currentDensity = float(obj.CurrentDensity_re_2.getValueAs("A/m^2")) + self.write.bodyForce(name, "Current Density 2", round(currentDensity, 6)) + if not obj.CurrentDensity_re_3_Disabled: + currentDensity = float(obj.CurrentDensity_re_3.getValueAs("A/m^2")) + self.write.bodyForce(name, "Current Density 3", round(currentDensity, 6)) + # imaginaries are only needed for harmonic equation + if equation.IsHarmonic: + if not obj.CurrentDensity_im_1_Disabled: + currentDensity = float(obj.CurrentDensity_im_1.getValueAs("A/m^2")) + self.write.bodyForce(name, "Current Density Im 1", round(currentDensity, 6)) + if not obj.CurrentDensity_im_2_Disabled: + currentDensity = float(obj.CurrentDensity_im_2.getValueAs("A/m^2")) + self.write.bodyForce(name, "Current Density Im 2", round(currentDensity, 6)) + if not obj.CurrentDensity_im_3_Disabled: + currentDensity = float(obj.CurrentDensity_im_3.getValueAs("A/m^2")) + self.write.bodyForce(name, "Current Density Im 3", round(currentDensity, 6)) + + if hasattr(obj, "Magnetization_im_1"): + # output only if magnetization is enabled and needed + if not obj.Magnetization_re_1_Disabled: + magnetization = float(obj.Magnetization_re_1.getValueAs("A/m")) + self.write.bodyForce(name, "Magnetization 1", magnetization) + if not obj.Magnetization_re_2_Disabled: + magnetization = float(obj.Magnetization_re_2.getValueAs("A/m")) + self.write.bodyForce(name, "Magnetization 2", magnetization) + if not obj.Magnetization_re_3_Disabled: + magnetization = float(obj.Magnetization_re_3.getValueAs("A/m")) + self.write.bodyForce(name, "Magnetization 3", magnetization) + # imaginaries are only needed for harmonic equation + if equation.IsHarmonic: + if not obj.Magnetization_im_1_Disabled: + magnetization = float(obj.Magnetization_im_1.getValueAs("A/m")) + self.write.bodyForce(name, "Magnetization Im 1", magnetization) + if not obj.Magnetization_im_2_Disabled: + magnetization = float(obj.Magnetization_im_2.getValueAs("A/m")) + self.write.bodyForce(name, "Magnetization Im 2", magnetization) + if not obj.Magnetization_im_3_Disabled: + magnetization = float(obj.Magnetization_im_3.getValueAs("A/m")) + self.write.bodyForce(name, "Magnetization Im 3", magnetization) + + if hasattr(obj, "PotentialEnabled"): + # check for PotentialEnabled not Potential since PotentialEnabled was + # added later and only with this the imaginary property is available + if obj.PotentialEnabled: + # output only if potential is enabled and needed + potential = float(obj.Potential.getValueAs("V")) + self.write.bodyForce(name, "Electric Potential", round(potential, 6)) + # imaginary is only needed for harmonic equation + if equation.IsHarmonic: + if not obj.AV_im_Disabled: + potential = float(obj.AV_im.getValueAs("V")) + self.write.bodyForce(name, "Electric Potential Im", round(potential, 6)) + + def handleMagnetodynamicBodyForces(self, bodies, equation): + # the current density can either be a body force or a boundary constraint + # therefore only output here if a solid is referenced + currentDensities = self.write.getMember("Fem::ConstraintCurrentDensity") + for obj in currentDensities: + if obj.References: + firstName = obj.References[0][1][0] + firstName = firstName.rstrip('0123456789') + if firstName == "Solid": + for name in obj.References[0][1]: + self._outputMagnetodynamicBodyForce(obj, name, equation) + self.write.handled(obj) + else: + # if there is only one current density without a reference, + # add it to all bodies + if len(currentDensities) == 1: + for name in bodies: + self._outputMagnetodynamicBodyForce(obj, name, equation) + else: + raise general_writer.WriteError( + "Several current density constraints found without reference to a body.\n" + "Please set a body for each current density constraint." + ) + self.write.handled(obj) + + magnetizations = self.write.getMember("Fem::ConstraintMagnetization") + for obj in magnetizations: + if obj.References: + for name in obj.References[0][1]: + self._outputMagnetodynamicBodyForce(obj, name, equation) + self.write.handled(obj) + else: + # if there is only one magnetization without a reference, + # add it to all bodies + if len(magnetizations) == 1: + for name in bodies: + self._outputMagnetodynamicBodyForce(obj, name, equation) + else: + raise general_writer.WriteError( + "Several magnetization constraints found without reference to a body.\n" + "Please set a body for each current density constraint." + ) + self.write.handled(obj) + + # the potential can either be a body force or a boundary constraint + # therefore only output here if a solid is referenced + potentials = self.write.getMember("Fem::ConstraintElectrostaticPotential") + for obj in potentials: + if obj.References: + firstName = obj.References[0][1][0] + firstName = firstName.rstrip('0123456789') + if firstName == "Solid": + for name in obj.References[0][1]: + # output only if potentiual is enabled and needed + self._outputMagnetodynamicBodyForce(obj, name, equation) + self.write.handled(obj) + + def _outputMagnetodynamicBndConditions(self, obj, name, equation): + if hasattr(obj, "CurrentDensity_re_1"): + # output only if current density is enabled and needed + if not obj.CurrentDensity_re_1_Disabled: + currentDensity = float(obj.CurrentDensity_re_1.getValueAs("A/m^2")) + self.write.boundary(name, "Current Density 1", round(currentDensity, 6)) + # imaginaries are only needed for harmonic equation + if equation.IsHarmonic: + if not obj.CurrentDensity_im_1_Disabled: + currentDensity = float(obj.CurrentDensity_im_1.getValueAs("A/m^2")) + self.write.boundary(name, "Current Density Im 1", round(currentDensity, 6)) + + if hasattr(obj, "PotentialEnabled"): + # check for PotentialEnabled not Potential since PotentialEnabled was + # added later and only with this the vectorial properties are available + if obj.PotentialEnabled: + potential = float(obj.Potential.getValueAs("V")) + if equation.IsHarmonic: + self.write.boundary(name, "AV re", round(potential, 6)) + else: + self.write.boundary(name, "AV", round(potential, 6)) + if not obj.AV_re_1_Disabled: + potential = float(obj.AV_re_1.getValueAs("V")) + if equation.IsHarmonic: + self.write.boundary(name, "AV re {e} 1", round(potential, 6)) + else: + self.write.boundary(name, "AV {e} 1", round(potential, 6)) + if not obj.AV_re_2_Disabled: + potential = float(obj.AV_re_2.getValueAs("V")) + if equation.IsHarmonic: + self.write.boundary(name, "AV re {e} 2", round(potential, 6)) + else: + self.write.boundary(name, "AV {e} 2", round(potential, 6)) + if not obj.AV_re_3_Disabled: + potential = float(obj.AV_re_3.getValueAs("V")) + if equation.IsHarmonic: + self.write.boundary(name, "AV re {e} 3", round(potential, 6)) + else: + self.write.boundary(name, "AV {e} 3", round(potential, 6)) + # imaginaries are only needed for harmonic equation + if equation.IsHarmonic: + if not obj.AV_im_Disabled: + potential = float(obj.AV_im.getValueAs("V")) + self.write.boundary(name, "AV im", round(potential, 6)) + if not obj.AV_im_1_Disabled: + potential = float(obj.AV_im_1.getValueAs("V")) + self.write.boundary(name, "AV im {e} 1", round(potential, 6)) + if not obj.AV_im_2_Disabled: + potential = float(obj.AV_im_2.getValueAs("V")) + self.write.boundary(name, "AV im {e} 2", round(potential, 6)) + if not obj.AV_im_3_Disabled: + potential = float(obj.AV_im_3.getValueAs("V")) + self.write.boundary(name, "AV im {e} 3", round(potential, 6)) + + def handleMagnetodynamicBndConditions(self, equation): + # the current density can either be a body force or a boundary constraint + # therefore only output here if a face is referenced + currentDensities = self.write.getMember("Fem::ConstraintCurrentDensity") + for obj in currentDensities: + if obj.References: + firstName = obj.References[0][1][0] + firstName = firstName.rstrip('0123456789') + if firstName == "Face": + for name in obj.References[0][1]: + self._outputMagnetodynamicBndConditions(obj, name, equation) + self.write.handled(obj) + + # the potential can either be a body force or a boundary constraint + # therefore only output here if a face is referenced + potentials = self.write.getMember("Fem::ConstraintElectrostaticPotential") + if len(potentials) == 0: + raise general_writer.WriteError( + "The Magnetodynamic equation needs at least one ElectrostaticPotential" + "constraint." + ) + for obj in potentials: + if obj.References: + firstName = obj.References[0][1][0] + firstName = firstName.rstrip('0123456789') + if firstName == "Face": + for name in obj.References[0][1]: + # output the FreeCAD label as comment + if obj.Label: + self.write.boundary(name, "! FreeCAD Name", obj.Label) + # output only if potentiual is enabled and needed + self._outputMagnetodynamicBndConditions(obj, name, equation) + self.write.handled(obj) + +## @} diff --git a/src/Mod/Fem/femsolver/elmer/writer.py b/src/Mod/Fem/femsolver/elmer/writer.py index 6872a91999..96f4cadab5 100644 --- a/src/Mod/Fem/femsolver/elmer/writer.py +++ b/src/Mod/Fem/femsolver/elmer/writer.py @@ -54,6 +54,7 @@ from .equations import electrostatic_writer as ES_writer from .equations import flow_writer from .equations import flux_writer from .equations import heat_writer +from .equations import magnetodynamic_writer as MgDyn_writer from .equations import magnetodynamic2D_writer as MgDyn2D_writer @@ -99,6 +100,7 @@ class Writer(object): self._handleHeat() self._handleFlow() self._handleFlux() + self._handleMagnetodynamic() self._handleMagnetodynamic2D() self._addOutputSolver() @@ -534,6 +536,30 @@ class Writer(object): HeatW.handleHeatBodyForces(activeIn) HeatW.handleHeatMaterial(activeIn) + #------------------------------------------------------------------------------------------- + # Magnetodynamic + + def _handleMagnetodynamic(self): + MgDyn = MgDyn_writer.MgDynwriter(self, self.solver) + activeIn = [] + for equation in self.solver.Group: + if femutils.is_of_type(equation, "Fem::EquationElmerMagnetodynamic"): + if equation.References: + activeIn = equation.References[0][1] + else: + activeIn = self.getAllBodies() + + solverSection = MgDyn.getMagnetodynamicSolver(equation) + solverPostSection = MgDyn.getMagnetodynamicSolverPost(equation) + for body in activeIn: + self._addSolver(body, solverSection) + self._addSolver(body, solverPostSection) + if activeIn: + MgDyn.handleMagnetodynamicConstants() + MgDyn.handleMagnetodynamicBndConditions(equation) + MgDyn.handleMagnetodynamicBodyForces(activeIn, equation) + MgDyn.handleMagnetodynamicMaterial(activeIn) + #------------------------------------------------------------------------------------------- # Magnetodynamic2D diff --git a/src/Mod/Fem/femsolver/equationbase.py b/src/Mod/Fem/femsolver/equationbase.py index e4fac3c488..a540c4165a 100644 --- a/src/Mod/Fem/femsolver/equationbase.py +++ b/src/Mod/Fem/femsolver/equationbase.py @@ -129,6 +129,16 @@ class HeatViewProxy(BaseViewProxy): return ":/icons/FEM_EquationHeat.svg" +class MagnetodynamicProxy(BaseProxy): + pass + + +class MagnetodynamicViewProxy(BaseViewProxy): + + def getIcon(self): + return ":/icons/FEM_EquationMagnetodynamic.svg" + + class Magnetodynamic2DProxy(BaseProxy): pass diff --git a/src/Mod/Fem/femtest/app/test_object.py b/src/Mod/Fem/femtest/app/test_object.py index a5a5d012c6..0f150527f1 100644 --- a/src/Mod/Fem/femtest/app/test_object.py +++ b/src/Mod/Fem/femtest/app/test_object.py @@ -84,14 +84,14 @@ class TestObjectCreate(unittest.TestCase): # thus they are not added to the analysis group ATM # https://forum.freecadweb.org/viewtopic.php?t=25283 # thus they should not be counted - # solver children: equations --> 7 + # solver children: equations --> 8 # gmsh mesh children: group, region, boundary layer --> 3 # resule children: mesh result --> 1 # post pipeline children: region, scalar, cut, wrap --> 4 # analysis itself is not in analysis group --> 1 - # thus: -16 + # thus: -17 - self.assertEqual(len(doc.Analysis.Group), count_defmake - 16) + self.assertEqual(len(doc.Analysis.Group), count_defmake - 17) self.assertEqual(len(doc.Objects), count_defmake) fcc_print("doc objects count: {}, method: {}".format( @@ -374,6 +374,10 @@ class TestObjectType(unittest.TestCase): "Fem::EquationElmerMagnetodynamic2D", type_of_obj(ObjectsFem.makeEquationMagnetodynamic2D(doc, solverelmer)) ) + self.assertEqual( + "Fem::EquationElmerMagnetodynamic", + type_of_obj(ObjectsFem.makeEquationMagnetodynamic(doc, solverelmer)) + ) fcc_print("doc objects count: {}, method: {}".format( len(doc.Objects), @@ -613,6 +617,10 @@ class TestObjectType(unittest.TestCase): ObjectsFem.makeEquationMagnetodynamic2D(doc, solverelmer), "Fem::EquationElmerMagnetodynamic2D" )) + self.assertTrue(is_of_type( + ObjectsFem.makeEquationMagnetodynamic(doc, solverelmer), + "Fem::EquationElmerMagnetodynamic" + )) fcc_print("doc objects count: {}, method: {}".format( len(doc.Objects), @@ -1453,7 +1461,7 @@ class TestObjectType(unittest.TestCase): "Fem::EquationElmerHeat" )) - # EquationElmerMagnetodynamic2D + # EquationElmerMagnetodynamic2D equation_magnetodynamic2D = ObjectsFem.makeEquationMagnetodynamic2D(doc, solver_elmer) self.assertTrue(is_derived_from( equation_magnetodynamic2D, @@ -1468,6 +1476,21 @@ class TestObjectType(unittest.TestCase): "Fem::EquationElmerMagnetodynamic2D" )) + # EquationElmerMagnetodynamic + equation_magnetodynamic = ObjectsFem.makeEquationMagnetodynamic(doc, solver_elmer) + self.assertTrue(is_derived_from( + equation_magnetodynamic, + "App::DocumentObject" + )) + self.assertTrue(is_derived_from( + equation_magnetodynamic, + "App::FeaturePython" + )) + self.assertTrue(is_derived_from( + equation_magnetodynamic, + "Fem::EquationElmerMagnetodynamic" + )) + fcc_print("doc objects count: {}, method: {}".format( len(doc.Objects), sys._getframe().f_code.co_name) @@ -1763,6 +1786,12 @@ class TestObjectType(unittest.TestCase): solverelmer ).isDerivedFrom("App::FeaturePython") ) + self.assertTrue( + ObjectsFem.makeEquationMagnetodynamic( + doc, + solverelmer + ).isDerivedFrom("App::FeaturePython") + ) fcc_print("doc objects count: {}, method: {}".format( len(doc.Objects), @@ -1845,6 +1874,7 @@ def create_all_fem_objects_doc( ObjectsFem.makeEquationFlux(doc, sol) ObjectsFem.makeEquationHeat(doc, sol) ObjectsFem.makeEquationMagnetodynamic2D(doc, sol) + ObjectsFem.makeEquationMagnetodynamic(doc, sol) doc.recompute() diff --git a/src/Mod/Fem/femtest/app/test_open.py b/src/Mod/Fem/femtest/app/test_open.py index a12b13c98b..9f952710a3 100644 --- a/src/Mod/Fem/femtest/app/test_open.py +++ b/src/Mod/Fem/femtest/app/test_open.py @@ -363,6 +363,11 @@ class TestObjectOpen(unittest.TestCase): type_of_obj(ObjectsFem.makeEquationMagnetodynamic2D(doc)) ) + self.assertEqual( + "Fem::EquationElmerMagnetodynamic", + type_of_obj(ObjectsFem.makeEquationMagnetodynamic(doc)) + ) + """ # code was generated by the following code from a document with all objects diff --git a/src/Mod/Fem/femtest/gui/test_open.py b/src/Mod/Fem/femtest/gui/test_open.py index 2ea21655c2..4f8db47664 100644 --- a/src/Mod/Fem/femtest/gui/test_open.py +++ b/src/Mod/Fem/femtest/gui/test_open.py @@ -338,6 +338,11 @@ class TestObjectOpen(unittest.TestCase): type_of_obj(ObjectsFem.makeEquationMagnetodynamic2D(doc)) ) + self.assertEqual( + "Fem::EquationElmerMagnetodynamic", + type_of_obj(ObjectsFem.makeEquationMagnetodynamic(doc)) + ) + """ # code was generated by the following code from a document with all objects