kindred/create
1
1
Fork 0
Code Issues 53 Pull Requests Actions Packages Projects 9 Releases 2 Wiki Activity

Fem: Add CalculiX writers for electrostatic constraints

Browse Source
This commit is contained in:
marioalexis
2025-05-07 18:49:05 -03:00
committed by Benjamin Nauck
parent eb82923869
commit 5e624a4775
9 changed files with 306 additions and 4 deletions
Download Patch File Download Diff File Expand all files Collapse all files

2
src/Mod/Fem/CMakeLists.txt
View File

@@ -248,6 +248,8 @@ SET(FemSolverCalculix_SRCS
femsolver/calculix/write_constraint_centrif.py
femsolver/calculix/write_constraint_contact.py
femsolver/calculix/write_constraint_displacement.py
femsolver/calculix/write_constraint_electricchargedensity.py
femsolver/calculix/write_constraint_electrostatic.py
femsolver/calculix/write_constraint_fixed.py
femsolver/calculix/write_constraint_fluidsection.py
femsolver/calculix/write_constraint_force.py

149
src/Mod/Fem/femsolver/calculix/write_constraint_electricchargedensity.py Normal file
View File

@@ -0,0 +1,149 @@
# SPDX-License-Identifier: LGPL-2.1-or-later
# ***************************************************************************
# * Copyright (c) 2025 Mario Passaglia <mpassaglia[at]cbc.uba.ar> *
# * *
# * This file is part of FreeCAD. *
# * *
# * FreeCAD is free software: you can redistribute it and/or modify it *
# * under the terms of the GNU Lesser General Public License as *
# * published by the Free Software Foundation, either version 2.1 of the *
# * License, or (at your option) any later version. *
# * *
# * FreeCAD 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 *
# * Lesser General Public License for more details. *
# * *
# * You should have received a copy of the GNU Lesser General Public *
# * License along with FreeCAD. If not, see *
# * <https://www.gnu.org/licenses/>. *
# * *
# ***************************************************************************
__title__ = "FreeCAD FEM calculix constraint electric charge density"
__author__ = "Mario Passaglia"
__url__ = "https://www.freecad.org"
def get_analysis_types():
return ["electromagnetic"]
def get_sets_name():
return "constraints_electricchargedensity_node_sets"
def get_constraint_title():
return "Electric charge density constraint applied"
def write_meshdata_constraint(f, femobj, den_obj, ccxwriter):
# print("MESHDATA: ", femobj)
if ccxwriter.solver_obj.ElectromagneticMode != "electrostatic":
return
if den_obj.Mode in ["Source", "Total Source"]:
msg, data = femobj["ChargeDensityElements"]
f.write(f"** {msg}\n")
f.write(f"*ELSET,ELSET={den_obj.Name}\n")
for ref, elem in data:
for e in elem:
f.write(f"{e},\n")
# if isinstance(femobj["FEMElements"], str):
# f.write("{}\n".format(femobj["FEMElements"]))
# else:
# for e in femobj["FEMElements"]:
# f.write(f"{e},\n")
# if femobj["Object"].BoundaryCondition == "Dirichlet":
# f.write(f"*NSET,NSET={den_obj.Name}\n")
# for n in femobj["Nodes"]:
# f.write(f"{n},\n")
def get_before_write_meshdata_constraint():
return ""
def get_after_write_meshdata_constraint():
return ""
def get_before_write_constraint():
return ""
def get_after_write_constraint():
return ""
def write_constraint(f, femobj, den_obj, ccxwriter):
# floats read from ccx should use {:.13G}, see comment in writer module
# if den_obj.BoundaryCondition == "Dirichlet":
# f.write("*BOUNDARY\n")
# f.write("{},11,11,{:.13G}\n".format(den_obj.Name, den_obj.Potential.getValueAs("mV").Value))
# f.write("\n")
if ccxwriter.solver_obj.ElectromagneticMode != "electrostatic":
return
match den_obj.Mode:
case "Source":
density = den_obj.SourceChargeDensity
case "Total Source":
density = den_obj.Proxy.get_total_source_density(den_obj)
case "Interface":
density = den_obj.InterfaceChargeDensity
case "Total Interface":
density = den_obj.Proxy.get_total_interface_density(den_obj)
if den_obj.Mode in ["Source", "Total Source"]:
f.write("*DFLUX\n")
f.write("{},BF,{:.13G}\n".format(den_obj.Name, density.getValueAs("C/mm^3").Value))
f.write("\n")
elif den_obj.Mode in ["Interface", "Total Interface"]:
# check internal interface
density = density.getValueAs("C/mm^2").Value
internal = _check_shared_interface(den_obj)
for feat, refs in femobj["ChargeDensityFaces"]:
f.write("** " + feat + "\n")
f.write("*DFLUX\n")
for ref in refs:
d = density
if ref[0] in internal:
d = density / 2
for face, fno in ref[1]:
f.write("{},S{},{:.13G}\n".format(face, fno, d))
f.write("\n")
def _check_shared_interface(den_obj):
"""
Check if reference is internal shared subshape
For example, shared face in compsolid
"""
internal = []
for o, sub in den_obj.References:
for elem in sub:
found = []
elem_i = o.getSubObject(elem)
if elem_i.ShapeType == "Face":
for s in o.Shape.Solids:
found.append(any([q.isSame(elem_i) for q in s.Faces]))
if sum(found) > 1:
internal.append((o, (elem,)))
if elem_i.ShapeType == "Edge":
for s in o.Shape.Faces:
found.append(any([q.isSame(elem_i) for q in s.Edges]))
if sum(found) > 1:
internal.append((o, (elem,)))
return internal

118
src/Mod/Fem/femsolver/calculix/write_constraint_electrostatic.py Normal file
View File

@@ -0,0 +1,118 @@
# SPDX-License-Identifier: LGPL-2.1-or-later
# ***************************************************************************
# * Copyright (c) 2025 Mario Passaglia <mpassaglia[at]cbc.uba.ar> *
# * *
# * This file is part of FreeCAD. *
# * *
# * FreeCAD is free software: you can redistribute it and/or modify it *
# * under the terms of the GNU Lesser General Public License as *
# * published by the Free Software Foundation, either version 2.1 of the *
# * License, or (at your option) any later version. *
# * *
# * FreeCAD 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 *
# * Lesser General Public License for more details. *
# * *
# * You should have received a copy of the GNU Lesser General Public *
# * License along with FreeCAD. If not, see *
# * <https://www.gnu.org/licenses/>. *
# * *
# ***************************************************************************
__title__ = "FreeCAD FEM calculix constraint electrostatic"
__author__ = "Mario Passaglia"
__url__ = "https://www.freecad.org"
import FreeCAD
def get_analysis_types():
return ["electromagnetic"]
def get_sets_name():
return "constraints_electrostaticpotential_node_sets"
def get_constraint_title():
return "Fixed electrostatic constraint applied"
def write_meshdata_constraint(f, femobj, pot_obj, ccxwriter):
if ccxwriter.solver_obj.ElectromagneticMode != "electrostatic":
return
if femobj["Object"].BoundaryCondition == "Dirichlet":
f.write(f"*NSET,NSET={pot_obj.Name}\n")
for n in femobj["Nodes"]:
f.write(f"{n},\n")
def get_before_write_meshdata_constraint():
return ""
def get_after_write_meshdata_constraint():
return ""
def get_before_write_constraint():
return ""
def get_after_write_constraint():
return ""
def write_constraint(f, femobj, pot_obj, ccxwriter):
if ccxwriter.solver_obj.ElectromagneticMode != "electrostatic":
return
# floats read from ccx should use {:.13G}, see comment in writer module
if pot_obj.BoundaryCondition == "Dirichlet":
f.write("*BOUNDARY\n")
f.write("{},11,11,{:.13G}\n".format(pot_obj.Name, pot_obj.Potential.getValueAs("mV").Value))
elif pot_obj.BoundaryCondition == "Neumann":
density = pot_obj.ElectricFluxDensity.getValueAs("C/mm^2").Value
# check internal interface
internal = _check_shared_interface(pot_obj)
for feat, refs in femobj["ElectricFluxFaces"]:
f.write("** " + feat + "\n")
f.write("*DFLUX\n")
for ref in refs:
d = density
if ref[0] in internal:
d = density / 2
for face, fno in ref[1]:
f.write("{},S{},{:.13G}\n".format(face, fno, d))
f.write("\n")
def _check_shared_interface(pot_obj):
"""
Check if reference is internal shared subshape
For example, shared face in compsolid
"""
internal = []
for o, sub in pot_obj.References:
for elem in sub:
found = []
elem_i = o.getSubObject(elem)
if elem_i.ShapeType == "Face":
for s in o.Shape.Solids:
found.append(any([q.isSame(elem_i) for q in s.Faces]))
if sum(found) > 1:
internal.append((o, (elem,)))
if elem_i.ShapeType == "Edge":
for s in o.Shape.Faces:
found.append(any([q.isSame(elem_i) for q in s.Edges]))
if sum(found) > 1:
internal.append((o, (elem,)))
return internal

15
src/Mod/Fem/femsolver/calculix/write_femelement_material.py
View File

@@ -27,6 +27,7 @@ __url__ = "https://www.freecad.org"
import FreeCAD
from femtools import constants
def write_femelement_material(f, ccxwriter):
@@ -92,7 +93,13 @@ def write_femelement_material(f, ccxwriter):
KV = FreeCAD.Units.Quantity(mat_obj.Material["KinematicViscosity"])
KV_in_mm2s = KV.getValueAs("mm^2/s").Value
DV_in_tmms = KV_in_mm2s * density_in_tonne_per_mm3
if (
ccxwriter.analysis_type == "electromagnetic"
and ccxwriter.solver_obj.ElectromagneticMode == "electrostatic"
):
rel_perm = FreeCAD.Units.Quantity(mat_obj.Material["RelativePermittivity"]).Value
vacuum_perm = FreeCAD.Units.Quantity(constants.vacuum_permittivity())
abs_perm = vacuum_perm.getValueAs("C/(mV*mm)").Value * rel_perm
# write material properties
f.write(f"** FreeCAD material name: {mat_info_name}\n")
f.write(f"** {mat_label}\n")
@@ -114,6 +121,12 @@ def write_femelement_material(f, ccxwriter):
elif mat_obj.Category == "Fluid":
f.write("*FLUID CONSTANTS\n")
f.write(f"{SH_in_JkgK:.13G},{DV_in_tmms:.13G}\n")
if (
ccxwriter.analysis_type == "electromagnetic"
and ccxwriter.solver_obj.ElectromagneticMode == "electrostatic"
):
f.write("*CONDUCTIVITY\n")
f.write(f"{abs_perm:.13G}\n")
# nonlinear material properties
if ccxwriter.solver_obj.MaterialNonlinearity == "nonlinear":

6
src/Mod/Fem/femsolver/calculix/write_step_equation.py
View File

@@ -78,6 +78,8 @@ def write_step_equation(f, ccxwriter):
analysis_type = "*NO ANALYSIS"
elif ccxwriter.analysis_type == "buckling":
analysis_type = "*BUCKLE"
elif ccxwriter.analysis_type == "electromagnetic":
analysis_type = "*HEAT TRANSFER, STEADY STATE"
# analysis line --> solver type
# https://forum.freecad.org/viewtopic.php?f=18&t=43178
if ccxwriter.solver_obj.MatrixSolverType == "default":
@@ -95,9 +97,7 @@ def write_step_equation(f, ccxwriter):
# analysis line --> user defined incrementations --> parameter DIRECT
# --> completely switch off ccx automatic incrementation
if ccxwriter.solver_obj.IterationsUserDefinedIncrementations:
if ccxwriter.analysis_type == "static":
analysis_type += ", DIRECT"
elif ccxwriter.analysis_type == "thermomech":
if ccxwriter.analysis_type in ["static", "thermomech", "electromagnetic"]:
analysis_type += ", DIRECT"
elif ccxwriter.analysis_type == "frequency":
FreeCAD.Console.PrintMessage(

4
src/Mod/Fem/femsolver/calculix/write_step_output.py
View File

@@ -47,6 +47,8 @@ def write_step_output(f, ccxwriter):
f.write("U, NT\n")
else:
f.write("MF, PS\n")
elif ccxwriter.analysis_type == "electromagnetic":
f.write("NT\n")
else:
f.write("U\n")
if not ccxwriter.member.geos_fluidsection:
@@ -56,6 +58,8 @@ def write_step_output(f, ccxwriter):
variables += ", HFL"
if ccxwriter.solver_obj.MaterialNonlinearity == "nonlinear":
variables += ", PEEQ"
if ccxwriter.analysis_type == "electromagnetic":
variables = "HFL"
f.write(variables + "\n")

10
src/Mod/Fem/femsolver/calculix/writer.py
View File

@@ -53,6 +53,8 @@ from . import write_constraint_selfweight as con_selfweight
from . import write_constraint_temperature as con_temperature
from . import write_constraint_tie as con_tie
from . import write_constraint_transform as con_transform
from . import write_constraint_electricchargedensity as con_electricchargedensity
from . import write_constraint_electrostatic as con_electrostatic
from . import write_femelement_geometry
from . import write_femelement_material
from . import write_femelement_matgeosets
@@ -157,6 +159,10 @@ class FemInputWriterCcx(writerbase.FemInputWriter):
self.write_constraints_meshsets(inpfile, self.member.cons_planerotation, con_planerotation)
self.write_constraints_meshsets(inpfile, self.member.cons_transform, con_transform)
self.write_constraints_meshsets(inpfile, self.member.cons_temperature, con_temperature)
self.write_constraints_meshsets(
inpfile, self.member.cons_electricchargedensity, con_electricchargedensity
)
self.write_constraints_meshsets(inpfile, self.member.cons_electrostatic, con_electrostatic)
# surface sets
self.write_constraints_meshsets(inpfile, self.member.cons_contact, con_contact)
@@ -194,6 +200,10 @@ class FemInputWriterCcx(writerbase.FemInputWriter):
self.write_constraints_meshsets(inpfile, self.member.cons_pressure, con_pressure)
self.write_constraints_propdata(inpfile, self.member.cons_temperature, con_temperature)
self.write_constraints_meshsets(inpfile, self.member.cons_heatflux, con_heatflux)
self.write_constraints_propdata(
inpfile, self.member.cons_electricchargedensity, con_electricchargedensity
)
self.write_constraints_propdata(inpfile, self.member.cons_electrostatic, con_electrostatic)
con_fluidsection.write_constraints_fluidsection(inpfile, self)
# output and step end

2
src/Mod/Fem/femsolver/writerbase.py
View File

@@ -128,6 +128,8 @@ class FemInputWriter:
self.temperature_objects = member.cons_temperature
self.tie_objects = member.cons_tie
self.transform_objects = member.cons_transform
self.electrostatic_objects = member.cons_electrostatic
self.electricchargedensity_objects = member.cons_electricchargedensity
# meshdatagetter, for compatibility, same with all getter methods
self.meshdatagetter = meshsetsgetter.MeshSetsGetter(

4
src/Mod/Fem/femtools/membertools.py
View File

@@ -282,6 +282,10 @@ class AnalysisMember:
self.cons_temperature = self.get_several_member("Fem::ConstraintTemperature")
self.cons_tie = self.get_several_member("Fem::ConstraintTie")
self.cons_transform = self.get_several_member("Fem::ConstraintTransform")
self.cons_electrostatic = self.get_several_member("Fem::ConstraintElectrostaticPotential")
self.cons_electricchargedensity = self.get_several_member(
"Fem::ConstraintElectricChargeDensity"
)
def get_several_member(self, t):
return get_several_member(self.analysis, t)