diff --git a/anisotropy/anisotropy.py b/anisotropy/anisotropy.py index f034287..6fa3508 100644 --- a/anisotropy/anisotropy.py +++ b/anisotropy/anisotropy.py @@ -2,12 +2,6 @@ import os, sys import time from datetime import timedelta, datetime import shutil - -#ROOT = "/".join(__file__.split("/")[:-2]) -#sys.path.append(os.path.abspath(ROOT)) - -#from anisotropy.utils import struct -import toml import logging __version__ = "1.1" @@ -135,12 +129,6 @@ class Anisotropy(object): return "\n".join([ f"{ k }: { v }" for k, v in versions.items() ]) - def setupDB(self): - self.db = db.init(self.env["db_path"]) - - if not os.path.exists(self.env["db_path"]): - self.db.create_tables([Structure, Mesh]) - def evalEnvParameters(self): """ 'Uncompress' and eval environment parameters """ from math import sqrt @@ -230,8 +218,19 @@ class Anisotropy(object): entry["geometry"]["theta"] == theta: return entry - # SELECT * FROM structure LEFT OUTER JOIN mesh ON mesh.structure_id = structure.id WHERE name = "faceCentered" AND direction = "[1, 1, 1]" AND theta = 0.12; - # Structure.select().join(Mesh, JOIN.LEFT_OUTER, on = (Mesh.structure_id == Structure.id)).where(Structure.name == "simple", Structure.direction == "[1, 0, 0]", Structure.theta == 0.13).dicts().get() + + def setupDB(self): + self.db = db.init(self.env["db_path"]) + + if not os.path.exists(self.env["db_path"]): + self.db.create_tables([ + Structure, + Mesh, + SubMesh, + MeshResult + ]) + + @timer def updateDB(self): for entry in self.params: @@ -252,6 +251,10 @@ class Anisotropy(object): m = deepcopy(entry["mesh"]) + sm = deepcopy(entry.get("submesh", {})) + + mr = deepcopy(entry.get("meshResult", {})) + if not query.exists(): with self.db.atomic(): stab = Structure.create(**s) @@ -259,6 +262,12 @@ class Anisotropy(object): m.update(structure_id = stab) mtab = Mesh.create(**m) + sm.update(mesh_id = mtab) + smtab = SubMesh.create(**sm) + + mr.update(mesh_id = mtab) + mrtab = MeshResult.create(**mr) + else: with self.db.atomic(): (Structure.update(**s) @@ -275,11 +284,20 @@ class Anisotropy(object): ) .execute()) + (SubMesh.update(**sm) + .where( + Submesh.mesh_id == None # TODO: ??? + ) + .execute()) + + # TODO: for MeshResult + @timer def updateFromDB(self): squery = Structure.select().order_by(Structure.id) mquery = Mesh.select().order_by(Mesh.structure_id) + for s, m in zip(squery.dicts(), mquery.dicts()): name = s.pop("name") path = s.pop("path") diff --git a/anisotropy/bodyCentered.py b/anisotropy/bodyCentered.py index 973e535..6da3645 100644 --- a/anisotropy/bodyCentered.py +++ b/anisotropy/bodyCentered.py @@ -1,318 +1,222 @@ -#import salome -#salome.salome_init() - -import GEOM -from salome.geom import geomBuilder -geompy = geomBuilder.New() - from math import pi, sqrt -def bodyCenteredCubic(theta = 0.01, fillet = False, direction = [1, 0, 0]): - - ### - # Parameters - ## - L = 1.0 - r0 = L * sqrt(3) / 4 +class simple(object): + def __init__(self, **kwargs): - radius = r0 / (1 - theta) - xn, yn, zn = 3, 3, 3 - - length = 2 * r0 - width = L / 2 - diag = L * sqrt(2) - height = L + self.direction = kwargs.get("direction", [1, 0, 0]) + self.theta = kwargs.get("theta", 0.01) + self.L = kwargs.get("L", 1) + self.r0 = kwargs.get("r0", self.L * sqrt(3) / 4) + self.radius = kwargs.get("radius", self.r0 / (1 - self.theta)) + self.filletsEnabled = kwargs.get("filletsEnabled", False) + self.fillets = kwargs.get("fillets", 0) - point = [] - xl = sqrt(diag ** 2 + diag ** 2) * 0.5 - yw = xl - zh = height - C1, C2 = 0.3, 0.2 - theta1, theta2 = 0.01, 0.18 - Cf = C1 + (C2 - C1) / (theta2 - theta1) * (theta - theta1) - delta = 0.02 - filletradius = delta - Cf * (radius - r0) - - scale = 100 - oo = geompy.MakeVertex(0, 0, 0) - spos1 = (0, 0, 0) - spos2 = (0, 0, 0) + def build(self): + from salomepl import getGeom - ### - # Bounding box - ## - if direction == [1, 0, 0]: - sk = geompy.Sketcher3D() - sk.addPointsAbsolute(xl, 0, 0) - sk.addPointsAbsolute(0, yw, 0) - sk.addPointsAbsolute(0, yw, zh) - sk.addPointsAbsolute(xl, 0, zh) - sk.addPointsAbsolute(xl, 0, 0) + geompy = getGeom() - inletface = geompy.MakeFaceWires([sk.wire()], True) - vecflow = geompy.GetNormal(inletface) - cubic = geompy.MakePrismVecH(inletface, vecflow, diag) + ### + # Pore Cell + ## + if self.direction in [[1, 0, 0], [0, 0, 1]]: + ### + # Parameters + ## + xn, yn, zn = 3, 3, 3 + + length = 2 * self.r0 + width = self.L / 2 + diag = self.L * sqrt(2) + height = self.L - elif direction == [0, 0, 1]: - sk = geompy.Sketcher3D() - sk.addPointsAbsolute(0, yw, 0) - sk.addPointsAbsolute(xl, 0, 0) - sk.addPointsAbsolute(2 * xl, yw, 0) - sk.addPointsAbsolute(xl, 2 * yw, 0) - sk.addPointsAbsolute(0, yw, 0) + point = [] + xl = sqrt(diag ** 2 + diag ** 2) * 0.5 + yw = xl + zh = height - inletface = geompy.MakeFaceWires([sk.wire()], True) - vecflow = geompy.GetNormal(inletface) - cubic = geompy.MakePrismVecH(inletface, vecflow, zh) + scale = 100 + oo = geompy.MakeVertex(0, 0, 0) + spos1 = (0, 0, 0) + spos2 = (0, 0, 0) - - inletface = geompy.MakeScaleTransform(inletface, oo, scale) - cubic = geompy.MakeScaleTransform(cubic, oo, scale) + ### + # Bounding box + ## + if direction == [1, 0, 0]: + sk = geompy.Sketcher3D() + sk.addPointsAbsolute(xl, 0, 0) + sk.addPointsAbsolute(0, yw, 0) + sk.addPointsAbsolute(0, yw, zh) + sk.addPointsAbsolute(xl, 0, zh) + sk.addPointsAbsolute(xl, 0, 0) - faces = geompy.ExtractShapes(cubic, geompy.ShapeType["FACE"], False) - symetryface = [] + inletface = geompy.MakeFaceWires([sk.wire()], True) + vecflow = geompy.GetNormal(inletface) + poreCell = geompy.MakePrismVecH(inletface, vecflow, diag) - for face in faces: - norm = geompy.GetNormal(face) - angle = round(geompy.GetAngle(norm, vecflow), 0) + elif direction == [0, 0, 1]: + sk = geompy.Sketcher3D() + sk.addPointsAbsolute(0, yw, 0) + sk.addPointsAbsolute(xl, 0, 0) + sk.addPointsAbsolute(2 * xl, yw, 0) + sk.addPointsAbsolute(xl, 2 * yw, 0) + sk.addPointsAbsolute(0, yw, 0) - if (angle == 0 or angle == 180) and not face == inletface: - outletface = face + inletface = geompy.MakeFaceWires([sk.wire()], True) + vecflow = geompy.GetNormal(inletface) + poreCell = geompy.MakePrismVecH(inletface, vecflow, zh) + + + inletface = geompy.MakeScaleTransform(inletface, oo, scale) + poreCell = geompy.MakeScaleTransform(poreCell, oo, scale) + + faces = geompy.ExtractShapes(poreCell, geompy.ShapeType["FACE"], False) + symetryface = [] + + for face in faces: + norm = geompy.GetNormal(face) + angle = round(geompy.GetAngle(norm, vecflow), 0) + + if (angle == 0 or angle == 180) and not face == inletface: + outletface = face + else: + symetryface.append(face) + + elif self.direction == [1, 1, 1]: + ### + # Parameters + ## + xn, yn, zn = 3, 3, 3 + + length = 2 * self.r0 + width = self.L / 2 + diag = self.L * sqrt(2) + height = diag / 3 + + point = [] + xl, yw, zh = -self.L / 4, -self.L / 4, -self.L / 4 + point.append((self.L / 3 + xl, self.L / 3 + yw, 4 * self.L / 3 + zh)) + point.append((self.L + xl, 0 + yw, self.L + zh)) + point.append((4 * self.L / 3 + xl, self.L / 3 + yw, self.L / 3 + zh)) + point.append((self.L + xl, self.L + yw, 0 + zh)) + point.append((self.L / 3 + xl, 4 * self.L / 3 + yw, self.L / 3 + zh)) + point.append((0 + xl, self.L + yw, self.L + zh)) + point.append((self.L / 3 + xl, self.L / 3 + yw, 4 * self.L / 3 + zh)) + + scale = 100 + oo = geompy.MakeVertex(0, 0, 0) + spos1 = (0, 0, 0) + spos2 = (0, 0, 0) + + ### + # Bounding box + ## + sk = geompy.Sketcher3D() + + for p in point: + sk.addPointsAbsolute(*p) + + inletface = geompy.MakeFaceWires([sk.wire()], False) + vecflow = geompy.GetNormal(inletface) + poreCell = geompy.MakePrismVecH(inletface, vecflow, self.L * sqrt(3)) + + inletface = geompy.MakeScaleTransform(inletface, oo, scale) + poreCell = geompy.MakeScaleTransform(poreCell, oo, scale) + + faces = geompy.ExtractShapes(poreCell, geompy.ShapeType["FACE"], False) + symetryface = [] + + for face in faces: + norm = geompy.GetNormal(face) + angle = round(geompy.GetAngle(norm, vecflow), 0) + + if (angle == 0 or angle == 180) and not face == inletface: + outletface = face + + else: + symetryface.append(face) + else: - symetryface.append(face) - - ### - # Grains - ## - ox = geompy.MakeVectorDXDYDZ(1, 0, 0) - oy = geompy.MakeVectorDXDYDZ(0, 1, 0) - oz = geompy.MakeVectorDXDYDZ(0, 0, 1) - xy = geompy.MakeVectorDXDYDZ(1, 1, 0) - xmy = geompy.MakeVectorDXDYDZ(1, -1, 0) - - grain = geompy.MakeSpherePntR(geompy.MakeVertex(*spos1), radius) - lattice1 = geompy.MakeMultiTranslation2D(grain, ox, L, xn, oy, L, yn) - lattice1 = geompy.MakeMultiTranslation1D(lattice1, oz, L, zn) - - #grain = geompy.MakeSpherePntR(geompy.MakeVertex(*spos2), radius) - #lattice2 = geompy.MakeMultiTranslation2D(grain, xy, length, xn + 1, xmy, length, yn + 1) - #lattice2 = geompy.MakeMultiTranslation1D(lattice2, oz, L, zn) - lattice2 = geompy.MakeTranslation(lattice1, 0.5 * L, 0.5 * L, 0.5 * L) - - grains = geompy.ExtractShapes(lattice1, geompy.ShapeType["SOLID"], True) - grains += geompy.ExtractShapes(lattice2, geompy.ShapeType["SOLID"], True) - grains = geompy.MakeFuseList(grains, False, False) - - grains = geompy.MakeScaleTransform(grains, oo, scale) - grainsOrigin = None - - if fillet: - grainsOrigin = geompy.MakeScaleTransform(grains, oo, 1 / scale) - grains = geompy.MakeFilletAll(grains, filletradius * scale) - - ### - # Groups - ## - shape = geompy.MakeCutList(cubic, [grains]) - shape = geompy.MakeScaleTransform(shape, oo, 1 / scale, theName = "bodyCenteredCubic") - - sall = geompy.CreateGroup(shape, geompy.ShapeType["FACE"]) - geompy.UnionIDs(sall, - geompy.SubShapeAllIDs(shape, geompy.ShapeType["FACE"])) - - inlet = geompy.CreateGroup(shape, geompy.ShapeType["FACE"], theName = "inlet") - inletshape = geompy.MakeCutList(inletface, [grains]) - inletshape = geompy.MakeScaleTransform(inletshape, oo, 1 / scale) - geompy.UnionList(inlet, geompy.SubShapeAll( - geompy.GetInPlace(shape, inletshape, True), geompy.ShapeType["FACE"])) - - outlet = geompy.CreateGroup(shape, geompy.ShapeType["FACE"], theName = "outlet") - outletshape = geompy.MakeCutList(outletface, [grains]) - outletshape = geompy.MakeScaleTransform(outletshape, oo, 1 / scale) - geompy.UnionList(outlet, geompy.SubShapeAll( - geompy.GetInPlace(shape, outletshape, True), geompy.ShapeType["FACE"])) - - symetry = [] - for (n, face) in enumerate(symetryface): - name = "symetry" + str(n) - symetry.append(geompy.CreateGroup(shape, geompy.ShapeType["FACE"], theName = name)) - symetryshape = geompy.MakeCutList(face, [grains]) - symetryshape = geompy.MakeScaleTransform(symetryshape, oo, 1 / scale) - geompy.UnionList(symetry[n], geompy.SubShapeAll( - geompy.GetInPlace(shape, symetryshape, True), geompy.ShapeType["FACE"])) - - groups = [] - groups.append(inlet) - groups.append(outlet) - groups.extend(symetry) - - if fillet: - strips = geompy.CreateGroup(shape, geompy.ShapeType["FACE"], theName = "strips") - shapeShell = geompy.ExtractShapes(shape, geompy.ShapeType["SHELL"], True) - stripsShape = geompy.MakeCutList(shapeShell[0], groups + [grainsOrigin]) - geompy.UnionList(strips, geompy.SubShapeAll( - geompy.GetInPlace(shape, stripsShape, True), geompy.ShapeType["FACE"])) - groups.append(strips) - - wall = geompy.CutListOfGroups([sall], groups, theName = "wall") - groups.append(wall) - - return shape, groups - - -def bodyCenteredHexagonalPrism(theta = 0.01, fillet = False, direction = [1, 1, 1]): - - ### - # Parameters - ## - L = 1.0 - r0 = L * sqrt(3) / 4 - - radius = r0 / (1 - theta) - xn, yn, zn = 3, 3, 3 - - length = 2 * r0 - width = L / 2 - diag = L * sqrt(2) - height = diag / 3 - - point = [] - xl, yw, zh = -L / 4, -L / 4, -L / 4 - point.append((L / 3 + xl, L / 3 + yw, 4 * L / 3 + zh)) - point.append((L + xl, 0 + yw, L + zh)) - point.append((4 * L / 3 + xl, L / 3 + yw, L / 3 + zh)) - point.append((L + xl, L + yw, 0 + zh)) - point.append((L / 3 + xl, 4 * L / 3 + yw, L / 3 + zh)) - point.append((0 + xl, L + yw, L + zh)) - point.append((L / 3 + xl, L / 3 + yw, 4 * L / 3 + zh)) - - C1, C2 = 0.3, 0.2 - theta1, theta2 = 0.01, 0.18 - Cf = C1 + (C2 - C1) / (theta2 - theta1) * (theta - theta1) - delta = 0.02 - filletradius = delta - Cf * (radius - r0) - - scale = 100 - oo = geompy.MakeVertex(0, 0, 0) - spos1 = (0, 0, 0) - spos2 = (0, 0, 0) - - ### - # Bounding box - ## - sk = geompy.Sketcher3D() - - for p in point: - sk.addPointsAbsolute(*p) - - inletface = geompy.MakeFaceWires([sk.wire()], False) - vecflow = geompy.GetNormal(inletface) - hexagonPrism = geompy.MakePrismVecH(inletface, vecflow, L * sqrt(3)) - - inletface = geompy.MakeScaleTransform(inletface, oo, scale) - hexagonPrism = geompy.MakeScaleTransform(hexagonPrism, oo, scale) - - faces = geompy.ExtractShapes(hexagonPrism, geompy.ShapeType["FACE"], False) - symetryface = [] - - for face in faces: - norm = geompy.GetNormal(face) - angle = round(geompy.GetAngle(norm, vecflow), 0) - - if (angle == 0 or angle == 180) and not face == inletface: - outletface = face + raise Exception(f"Direction { self.direction } is not implemented") - else: - symetryface.append(face) - - ### - # Grains - ## - ox = geompy.MakeVectorDXDYDZ(1, 0, 0) - oy = geompy.MakeVectorDXDYDZ(0, 1, 0) - oz = geompy.MakeVectorDXDYDZ(0, 0, 1) - xy = geompy.MakeVectorDXDYDZ(1, 1, 0) - xmy = geompy.MakeVectorDXDYDZ(1, -1, 0) + ### + # Grains + ## + ox = geompy.MakeVectorDXDYDZ(1, 0, 0) + oy = geompy.MakeVectorDXDYDZ(0, 1, 0) + oz = geompy.MakeVectorDXDYDZ(0, 0, 1) + xy = geompy.MakeVectorDXDYDZ(1, 1, 0) + xmy = geompy.MakeVectorDXDYDZ(1, -1, 0) - grain = geompy.MakeSpherePntR(geompy.MakeVertex(*spos1), radius) - lattice1 = geompy.MakeMultiTranslation2D(grain, ox, L, xn, oy, L, yn) - lattice1 = geompy.MakeMultiTranslation1D(lattice1, oz, L, zn) + grain = geompy.MakeSpherePntR(geompy.MakeVertex(*spos1), self.radius) + lattice1 = geompy.MakeMultiTranslation2D(grain, ox, self.L, xn, oy, self.L, yn) + lattice1 = geompy.MakeMultiTranslation1D(lattice1, oz, self.L, zn) - #grain = geompy.MakeSpherePntR(geompy.MakeVertex(*spos2), radius) - #lattice2 = geompy.MakeMultiTranslation2D(grain, xy, length, xn + 1, xmy, length, yn + 1) - #lattice2 = geompy.MakeMultiTranslation1D(lattice2, oz, L, zn) - lattice2 = geompy.MakeTranslation(lattice1, 0.5 * L, 0.5 * L, 0.5 * L) - - grains = geompy.ExtractShapes(lattice1, geompy.ShapeType["SOLID"], True) - grains += geompy.ExtractShapes(lattice2, geompy.ShapeType["SOLID"], True) - grains = geompy.MakeFuseList(grains, False, False) + #grain = geompy.MakeSpherePntR(geompy.MakeVertex(*spos2), radius) + #lattice2 = geompy.MakeMultiTranslation2D(grain, xy, length, xn + 1, xmy, length, yn + 1) + #lattice2 = geompy.MakeMultiTranslation1D(lattice2, oz, L, zn) + lattice2 = geompy.MakeTranslation(lattice1, 0.5 * self.L, 0.5 * self.L, 0.5 * self.L) + + grains = geompy.ExtractShapes(lattice1, geompy.ShapeType["SOLID"], True) + grains += geompy.ExtractShapes(lattice2, geompy.ShapeType["SOLID"], True) + grains = geompy.MakeFuseList(grains, False, False) - grains = geompy.MakeScaleTransform(grains, oo, scale) - grainsOrigin = None + grains = geompy.MakeScaleTransform(grains, oo, scale) + grainsOrigin = None - if fillet: - grainsOrigin = geompy.MakeScaleTransform(grains, oo, 1 / scale) - grains = geompy.MakeFilletAll(grains, filletradius * scale) + if self.filletsEnabled: + grainsOrigin = geompy.MakeScaleTransform(grains, oo, 1 / scale) + grains = geompy.MakeFilletAll(grains, self.fillets * scale) - ### - # Groups - ## - shape = geompy.MakeCutList(hexagonPrism, [grains]) - shape = geompy.MakeScaleTransform(shape, oo, 1 / scale, theName = "bodyCenteredCubic") - - sall = geompy.CreateGroup(shape, geompy.ShapeType["FACE"]) - geompy.UnionIDs(sall, - geompy.SubShapeAllIDs(shape, geompy.ShapeType["FACE"])) + ### + # Groups + ## + shape = geompy.MakeCutList(poreCell, [grains]) + shape = geompy.MakeScaleTransform(shape, oo, 1 / scale, theName = "bodyCenteredCubic") - inlet = geompy.CreateGroup(shape, geompy.ShapeType["FACE"], theName = "inlet") - inletshape = geompy.MakeCutList(inletface, [grains]) - inletshape = geompy.MakeScaleTransform(inletshape, oo, 1 / scale) - geompy.UnionList(inlet, geompy.SubShapeAll( - geompy.GetInPlace(shape, inletshape, True), geompy.ShapeType["FACE"])) + sall = geompy.CreateGroup(shape, geompy.ShapeType["FACE"]) + geompy.UnionIDs(sall, + geompy.SubShapeAllIDs(shape, geompy.ShapeType["FACE"])) - outlet = geompy.CreateGroup(shape, geompy.ShapeType["FACE"], theName = "outlet") - outletshape = geompy.MakeCutList(outletface, [grains]) - outletshape = geompy.MakeScaleTransform(outletshape, oo, 1 / scale) - geompy.UnionList(outlet, geompy.SubShapeAll( - geompy.GetInPlace(shape, outletshape, True), geompy.ShapeType["FACE"])) - - symetry = [] - for (n, face) in enumerate(symetryface): - name = "symetry" + str(n) - symetry.append(geompy.CreateGroup(shape, geompy.ShapeType["FACE"], theName = name)) - symetryshape = geompy.MakeCutList(face, [grains]) - symetryshape = geompy.MakeScaleTransform(symetryshape, oo, 1 / scale) - geompy.UnionList(symetry[n], geompy.SubShapeAll( - geompy.GetInPlace(shape, symetryshape, True), geompy.ShapeType["FACE"])) + inlet = geompy.CreateGroup(shape, geompy.ShapeType["FACE"], theName = "inlet") + inletshape = geompy.MakeCutList(inletface, [grains]) + inletshape = geompy.MakeScaleTransform(inletshape, oo, 1 / scale) + geompy.UnionList(inlet, geompy.SubShapeAll( + geompy.GetInPlace(shape, inletshape, True), geompy.ShapeType["FACE"])) - groups = [] - groups.append(inlet) - groups.append(outlet) - groups.extend(symetry) - - if fillet: - strips = geompy.CreateGroup(shape, geompy.ShapeType["FACE"], theName = "strips") - shapeShell = geompy.ExtractShapes(shape, geompy.ShapeType["SHELL"], True) - stripsShape = geompy.MakeCutList(shapeShell[0], groups + [grainsOrigin]) - geompy.UnionList(strips, geompy.SubShapeAll( - geompy.GetInPlace(shape, stripsShape, True), geompy.ShapeType["FACE"])) - groups.append(strips) - - wall = geompy.CutListOfGroups([sall], groups, theName = "wall") - groups.append(wall) + outlet = geompy.CreateGroup(shape, geompy.ShapeType["FACE"], theName = "outlet") + outletshape = geompy.MakeCutList(outletface, [grains]) + outletshape = geompy.MakeScaleTransform(outletshape, oo, 1 / scale) + geompy.UnionList(outlet, geompy.SubShapeAll( + geompy.GetInPlace(shape, outletshape, True), geompy.ShapeType["FACE"])) + + symetry = [] + for (n, face) in enumerate(symetryface): + name = "symetry" + str(n) + symetry.append(geompy.CreateGroup(shape, geompy.ShapeType["FACE"], theName = name)) + symetryshape = geompy.MakeCutList(face, [grains]) + symetryshape = geompy.MakeScaleTransform(symetryshape, oo, 1 / scale) + geompy.UnionList(symetry[n], geompy.SubShapeAll( + geompy.GetInPlace(shape, symetryshape, True), geompy.ShapeType["FACE"])) - return shape, groups + groups = [] + groups.append(inlet) + groups.append(outlet) + groups.extend(symetry) + if self.filletsEnabled: + strips = geompy.CreateGroup(shape, geompy.ShapeType["FACE"], theName = "strips") + shapeShell = geompy.ExtractShapes(shape, geompy.ShapeType["SHELL"], True) + stripsShape = geompy.MakeCutList(shapeShell[0], groups + [grainsOrigin]) + geompy.UnionList(strips, geompy.SubShapeAll( + geompy.GetInPlace(shape, stripsShape, True), geompy.ShapeType["FACE"])) + groups.append(strips) + + wall = geompy.CutListOfGroups([sall], groups, theName = "wall") + groups.append(wall) -def bodyCentered(theta, fillet, direction): - if direction in [[1, 0, 0], [0, 0, 1]]: - return bodyCenteredCubic(theta, fillet, direction) - - elif direction == [1, 1, 1]: - return bodyCenteredHexagonalPrism(theta, fillet, direction) - - else: - raise Exception("This direction is not implemented") + return shape, groups diff --git a/anisotropy/faceCentered.py b/anisotropy/faceCentered.py index 6952ce0..4c9214c 100644 --- a/anisotropy/faceCentered.py +++ b/anisotropy/faceCentered.py @@ -1,315 +1,222 @@ -#import salome -#salome.salome_init() - -import GEOM -from salome.geom import geomBuilder -geompy = geomBuilder.New() - from math import pi, sqrt -def faceCenteredCubic(theta = 0.01, fillet = False, direction = [1, 0, 0]): +class faceCentered(object): + def __init__(self, **kwargs): + + self.direction = kwargs.get("direction", [1, 0, 0]) + self.theta = kwargs.get("theta", 0.01) + self.L = kwargs.get("L", 1) + self.r0 = kwargs.get("r0", self.L * sqrt(2) / 4) + self.radius = kwargs.get("radius", self.r0 / (1 - self.theta)) + self.filletsEnabled = kwargs.get("filletsEnabled", False) + self.fillets = kwargs.get("fillets", 0) + + + def build(self): + from salomepl import getGeom + + geompy = getGeom() + + ### + # Pore Cell + ## + if self.direction in [[1, 0, 0], [0, 0, 1]]: + ### + # Parameters + ## + xn, yn, zn = 3, 3, 3 + + length = 2 * self.r0 + width = self.L / 2 + diag = self.L * sqrt(3) + height = diag / 3 + + point = [] + xl = sqrt(length ** 2 + length ** 2) * 0.5 + yw = xl + zh = width + + scale = 100 + oo = geompy.MakeVertex(0, 0, 0) + spos1 = (-width * (xn - 1), 0, -width * (zn - 2)) + spos2 = (-width * xn, 0, -width * (zn - 1)) + + ### + # Bounding box + ## + if self.direction == [1, 0, 0]: + sk = geompy.Sketcher3D() + sk.addPointsAbsolute(0, 0, -zh) + sk.addPointsAbsolute(-xl, yw, -zh) + sk.addPointsAbsolute(-xl, yw, zh) + sk.addPointsAbsolute(0, 0, zh) + sk.addPointsAbsolute(0, 0, -zh) + + inletface = geompy.MakeFaceWires([sk.wire()], True) + vecflow = geompy.GetNormal(inletface) + poreCell = geompy.MakePrismVecH(inletface, vecflow, length) + + elif self.direction == [0, 0, 1]: + sk = geompy.Sketcher3D() + sk.addPointsAbsolute(0, 0, -zh) + sk.addPointsAbsolute(xl, yw, -zh) + sk.addPointsAbsolute(0, 2 * yw, -zh) + sk.addPointsAbsolute(-xl, yw, -zh) + sk.addPointsAbsolute(0, 0, -zh) + + inletface = geompy.MakeFaceWires([sk.wire()], True) + vecflow = geompy.GetNormal(inletface) + poreCell = geompy.MakePrismVecH(inletface, vecflow, 2 * zh) + + + inletface = geompy.MakeScaleTransform(inletface, oo, scale) + poreCell = geompy.MakeScaleTransform(poreCell, oo, scale) + + faces = geompy.ExtractShapes(poreCell, geompy.ShapeType["FACE"], False) + symetryface = [] + + for face in faces: + norm = geompy.GetNormal(face) + angle = round(geompy.GetAngle(norm, vecflow), 0) + + if (angle == 0 or angle == 180) and not face == inletface: + outletface = face + + else: + symetryface.append(face) - ### - # Parameters - ## - L = 1.0 - r0 = L * sqrt(2) / 4 + elif self.direction == [1, 1, 1]: + ### + # Parameters + ## + xn, yn, zn = 3, 3, 3 + + length = 2 * self.r0 + width = self.L / 2 + diag = self.L * sqrt(3) + height = diag / 3 - radius = r0 / (1 - theta) - xn, yn, zn = 3, 3, 3 - - length = 2 * r0 - width = L / 2 - diag = L * sqrt(3) - height = diag / 3 + point = [] + xl, yw, zh = -(xn - 2) * self.L / 3, -(yn - 2) * self.L / 3, -(zn - 2) * self.L / 3 + point.append((-2 * width / 3 + xl, -2 * width / 3 + yw, width / 3 + zh)) + point.append((0 + xl, -width + yw, 0 + zh)) + point.append((width / 3 + xl, -2 * width / 3 + yw, -2 * width / 3 + zh)) + point.append((0 + xl, 0 + yw, -width + zh)) + point.append((-2 * width / 3 + xl, width / 3 + yw, -2 * width / 3 + zh)) + point.append((-width + xl, 0 + yw, 0 + zh)) + point.append((-2 * width / 3 + xl, -2 * width / 3 + yw, width / 3 + zh)) - point = [] - xl = sqrt(length ** 2 + length ** 2) * 0.5 - yw = xl - zh = width + scale = 100 + oo = geompy.MakeVertex(0, 0, 0) + spos1 = (-width * (xn - 1), 0, -width * (zn - 2)) + spos2 = (-width * xn, 0, -width * (zn - 1)) - C1, C2 = 0.3, 0.2 - theta1, theta2 = 0.01, 0.13 - Cf = C1 + (C2 - C1) / (theta2 - theta1) * (theta - theta1) - delta = 0.012 - filletradius = delta - Cf * (radius - r0) - - scale = 100 - oo = geompy.MakeVertex(0, 0, 0) - spos1 = (-width * (xn - 1), 0, -width * (zn - 2)) - spos2 = (-width * xn, 0, -width * (zn - 1)) + ### + # Bounding box + ## + sk = geompy.Sketcher3D() - ### - # Bounding box - ## - if direction == [1, 0, 0]: - sk = geompy.Sketcher3D() - sk.addPointsAbsolute(0, 0, -zh) - sk.addPointsAbsolute(-xl, yw, -zh) - sk.addPointsAbsolute(-xl, yw, zh) - sk.addPointsAbsolute(0, 0, zh) - sk.addPointsAbsolute(0, 0, -zh) + for p in point: + sk.addPointsAbsolute(*p) + + inletface = geompy.MakeFaceWires([sk.wire()], False) + vecflow = geompy.GetNormal(inletface) + poreCell = geompy.MakePrismVecH(inletface, vecflow, self.L * sqrt(3)) - inletface = geompy.MakeFaceWires([sk.wire()], True) - vecflow = geompy.GetNormal(inletface) - cubic = geompy.MakePrismVecH(inletface, vecflow, length) + inletface = geompy.MakeScaleTransform(inletface, oo, scale) + poreCell = geompy.MakeScaleTransform(poreCell, oo, scale) - elif direction == [0, 0, 1]: - sk = geompy.Sketcher3D() - sk.addPointsAbsolute(0, 0, -zh) - sk.addPointsAbsolute(xl, yw, -zh) - sk.addPointsAbsolute(0, 2 * yw, -zh) - sk.addPointsAbsolute(-xl, yw, -zh) - sk.addPointsAbsolute(0, 0, -zh) + faces = geompy.ExtractShapes(poreCell, geompy.ShapeType["FACE"], False) + symetryface = [] - inletface = geompy.MakeFaceWires([sk.wire()], True) - vecflow = geompy.GetNormal(inletface) - cubic = geompy.MakePrismVecH(inletface, vecflow, 2 * zh) + for face in faces: + norm = geompy.GetNormal(face) + angle = round(geompy.GetAngle(norm, vecflow), 0) - - inletface = geompy.MakeScaleTransform(inletface, oo, scale) - cubic = geompy.MakeScaleTransform(cubic, oo, scale) + if (angle == 0 or angle == 180) and not face == inletface: + outletface = face + + else: + symetryface.append(face) - faces = geompy.ExtractShapes(cubic, geompy.ShapeType["FACE"], False) - symetryface = [] - - for face in faces: - norm = geompy.GetNormal(face) - angle = round(geompy.GetAngle(norm, vecflow), 0) - - if (angle == 0 or angle == 180) and not face == inletface: - outletface = face - else: - symetryface.append(face) - - ### - # Grains - ## - ox = geompy.MakeVectorDXDYDZ(1, 0, 0) - oy = geompy.MakeVectorDXDYDZ(0, 1, 0) - oz = geompy.MakeVectorDXDYDZ(0, 0, 1) - xy = geompy.MakeVectorDXDYDZ(1, 1, 0) - xmy = geompy.MakeVectorDXDYDZ(1, -1, 0) - - grain = geompy.MakeSpherePntR(geompy.MakeVertex(*spos1), radius) - lattice1 = geompy.MakeMultiTranslation2D(grain, xy, length, xn, xmy, length, yn) - lattice1 = geompy.MakeMultiTranslation1D(lattice1, oz, L, zn - 1) - - grain = geompy.MakeSpherePntR(geompy.MakeVertex(*spos2), radius) - lattice2 = geompy.MakeMultiTranslation2D(grain, xy, length, xn + 1, xmy, length, yn + 1) - lattice2 = geompy.MakeMultiTranslation1D(lattice2, oz, L, zn) - - grains = geompy.ExtractShapes(lattice1, geompy.ShapeType["SOLID"], True) - grains += geompy.ExtractShapes(lattice2, geompy.ShapeType["SOLID"], True) - grains = geompy.MakeFuseList(grains, False, False) - - grains = geompy.MakeScaleTransform(grains, oo, scale) - grainsOrigin = None - - if fillet: - grainsOrigin = geompy.MakeScaleTransform(grains, oo, 1 / scale) - grains = geompy.MakeFilletAll(grains, filletradius * scale) - - ### - # Groups - ## - shape = geompy.MakeCutList(cubic, [grains]) - shape = geompy.MakeScaleTransform(shape, oo, 1 / scale, theName = "faceCenteredCubic") - - sall = geompy.CreateGroup(shape, geompy.ShapeType["FACE"]) - geompy.UnionIDs(sall, - geompy.SubShapeAllIDs(shape, geompy.ShapeType["FACE"])) - - inlet = geompy.CreateGroup(shape, geompy.ShapeType["FACE"], theName = "inlet") - inletshape = geompy.MakeCutList(inletface, [grains]) - inletshape = geompy.MakeScaleTransform(inletshape, oo, 1 / scale) - geompy.UnionList(inlet, geompy.SubShapeAll( - geompy.GetInPlace(shape, inletshape, True), geompy.ShapeType["FACE"])) - - outlet = geompy.CreateGroup(shape, geompy.ShapeType["FACE"], theName = "outlet") - outletshape = geompy.MakeCutList(outletface, [grains]) - outletshape = geompy.MakeScaleTransform(outletshape, oo, 1 / scale) - geompy.UnionList(outlet, geompy.SubShapeAll( - geompy.GetInPlace(shape, outletshape, True), geompy.ShapeType["FACE"])) - - symetry = [] - for (n, face) in enumerate(symetryface): - name = "symetry" + str(n) - symetry.append(geompy.CreateGroup(shape, geompy.ShapeType["FACE"], theName = name)) - symetryshape = geompy.MakeCutList(face, [grains]) - symetryshape = geompy.MakeScaleTransform(symetryshape, oo, 1 / scale) - geompy.UnionList(symetry[n], geompy.SubShapeAll( - geompy.GetInPlace(shape, symetryshape, True), geompy.ShapeType["FACE"])) - - groups = [] - groups.append(inlet) - groups.append(outlet) - groups.extend(symetry) - - if fillet: - strips = geompy.CreateGroup(shape, geompy.ShapeType["FACE"], theName = "strips") - shapeShell = geompy.ExtractShapes(shape, geompy.ShapeType["SHELL"], True) - stripsShape = geompy.MakeCutList(shapeShell[0], groups + [grainsOrigin]) - geompy.UnionList(strips, geompy.SubShapeAll( - geompy.GetInPlace(shape, stripsShape, True), geompy.ShapeType["FACE"])) - groups.append(strips) - - wall = geompy.CutListOfGroups([sall], groups, theName = "wall") - groups.append(wall) - - return shape, groups - -def faceCenteredHexagonalPrism(theta = 0.01, fillet = False, direction = [1, 1, 1]): - - ### - # Parameters - ## - L = 1.0 - r0 = L * sqrt(2) / 4 - - radius = r0 / (1 - theta) - xn, yn, zn = 3, 3, 3 - - length = 2 * r0 - width = L / 2 - diag = L * sqrt(3) - height = diag / 3 - - point = [] - xl, yw, zh = -(xn - 2) * L / 3, -(yn - 2) * L / 3, -(zn - 2) * L / 3 - point.append((-2 * width / 3 + xl, -2 * width / 3 + yw, width / 3 + zh)) - point.append((0 + xl, -width + yw, 0 + zh)) - point.append((width / 3 + xl, -2 * width / 3 + yw, -2 * width / 3 + zh)) - point.append((0 + xl, 0 + yw, -width + zh)) - point.append((-2 * width / 3 + xl, width / 3 + yw, -2 * width / 3 + zh)) - point.append((-width + xl, 0 + yw, 0 + zh)) - point.append((-2 * width / 3 + xl, -2 * width / 3 + yw, width / 3 + zh)) - - C1, C2 = 0.3, 0.2 - theta1, theta2 = 0.01, 0.13 - Cf = C1 + (C2 - C1) / (theta2 - theta1) * (theta - theta1) - delta = 0.012 - filletradius = delta - Cf * (radius - r0) - - scale = 100 - oo = geompy.MakeVertex(0, 0, 0) - spos1 = (-width * (xn - 1), 0, -width * (zn - 2)) - spos2 = (-width * xn, 0, -width * (zn - 1)) - - ### - # Bounding box - ## - sk = geompy.Sketcher3D() - - for p in point: - sk.addPointsAbsolute(*p) - - inletface = geompy.MakeFaceWires([sk.wire()], False) - vecflow = geompy.GetNormal(inletface) - hexagonPrism = geompy.MakePrismVecH(inletface, vecflow, L * sqrt(3)) - - inletface = geompy.MakeScaleTransform(inletface, oo, scale) - hexagonPrism = geompy.MakeScaleTransform(hexagonPrism, oo, scale) - - faces = geompy.ExtractShapes(hexagonPrism, geompy.ShapeType["FACE"], False) - symetryface = [] - - for face in faces: - norm = geompy.GetNormal(face) - angle = round(geompy.GetAngle(norm, vecflow), 0) - - if (angle == 0 or angle == 180) and not face == inletface: - outletface = face + raise Exception(f"Direction { self.direction } is not implemented") - else: - symetryface.append(face) - - ### - # Grains - ## - ox = geompy.MakeVectorDXDYDZ(1, 0, 0) - oy = geompy.MakeVectorDXDYDZ(0, 1, 0) - oz = geompy.MakeVectorDXDYDZ(0, 0, 1) - xy = geompy.MakeVectorDXDYDZ(1, 1, 0) - xmy = geompy.MakeVectorDXDYDZ(1, -1, 0) + ### + # Grains + ## + ox = geompy.MakeVectorDXDYDZ(1, 0, 0) + oy = geompy.MakeVectorDXDYDZ(0, 1, 0) + oz = geompy.MakeVectorDXDYDZ(0, 0, 1) + xy = geompy.MakeVectorDXDYDZ(1, 1, 0) + xmy = geompy.MakeVectorDXDYDZ(1, -1, 0) - grain = geompy.MakeSpherePntR(geompy.MakeVertex(*spos1), radius) - lattice1 = geompy.MakeMultiTranslation2D(grain, xy, length, xn, xmy, length, yn) - lattice1 = geompy.MakeMultiTranslation1D(lattice1, oz, L, zn - 1) + grain = geompy.MakeSpherePntR(geompy.MakeVertex(*spos1), radius) + lattice1 = geompy.MakeMultiTranslation2D(grain, xy, length, xn, xmy, length, yn) + lattice1 = geompy.MakeMultiTranslation1D(lattice1, oz, L, zn - 1) - grain = geompy.MakeSpherePntR(geompy.MakeVertex(*spos2), radius) - lattice2 = geompy.MakeMultiTranslation2D(grain, xy, length, xn + 1, xmy, length, yn + 1) - lattice2 = geompy.MakeMultiTranslation1D(lattice2, oz, L, zn) - - grains = geompy.ExtractShapes(lattice1, geompy.ShapeType["SOLID"], True) - grains += geompy.ExtractShapes(lattice2, geompy.ShapeType["SOLID"], True) - grains = geompy.MakeFuseList(grains, False, False) + grain = geompy.MakeSpherePntR(geompy.MakeVertex(*spos2), radius) + lattice2 = geompy.MakeMultiTranslation2D(grain, xy, length, xn + 1, xmy, length, yn + 1) + lattice2 = geompy.MakeMultiTranslation1D(lattice2, oz, L, zn) + + grains = geompy.ExtractShapes(lattice1, geompy.ShapeType["SOLID"], True) + grains += geompy.ExtractShapes(lattice2, geompy.ShapeType["SOLID"], True) + grains = geompy.MakeFuseList(grains, False, False) - grains = geompy.MakeScaleTransform(grains, oo, scale) - grainsOrigin = None + grains = geompy.MakeScaleTransform(grains, oo, scale) + grainsOrigin = None - if fillet: - grainsOrigin = geompy.MakeScaleTransform(grains, oo, 1 / scale) - grains = geompy.MakeFilletAll(grains, filletradius * scale) + if self.filletsEnabled: + grainsOrigin = geompy.MakeScaleTransform(grains, oo, 1 / scale) + grains = geompy.MakeFilletAll(grains, self.fillets * scale) - ### - # Groups - ## - shape = geompy.MakeCutList(hexagonPrism, [grains]) - shape = geompy.MakeScaleTransform(shape, oo, 1 / scale, theName = "faceCenteredCubic") - - sall = geompy.CreateGroup(shape, geompy.ShapeType["FACE"]) - geompy.UnionIDs(sall, - geompy.SubShapeAllIDs(shape, geompy.ShapeType["FACE"])) + ### + # Groups + ## + shape = geompy.MakeCutList(poreCell, [grains]) + shape = geompy.MakeScaleTransform(shape, oo, 1 / scale, theName = "faceCentered") + + sall = geompy.CreateGroup(shape, geompy.ShapeType["FACE"]) + geompy.UnionIDs(sall, + geompy.SubShapeAllIDs(shape, geompy.ShapeType["FACE"])) - inlet = geompy.CreateGroup(shape, geompy.ShapeType["FACE"], theName = "inlet") - inletshape = geompy.MakeCutList(inletface, [grains]) - inletshape = geompy.MakeScaleTransform(inletshape, oo, 1 / scale) - geompy.UnionList(inlet, geompy.SubShapeAll( - geompy.GetInPlace(shape, inletshape, True), geompy.ShapeType["FACE"])) + inlet = geompy.CreateGroup(shape, geompy.ShapeType["FACE"], theName = "inlet") + inletshape = geompy.MakeCutList(inletface, [grains]) + inletshape = geompy.MakeScaleTransform(inletshape, oo, 1 / scale) + geompy.UnionList(inlet, geompy.SubShapeAll( + geompy.GetInPlace(shape, inletshape, True), geompy.ShapeType["FACE"])) - outlet = geompy.CreateGroup(shape, geompy.ShapeType["FACE"], theName = "outlet") - outletshape = geompy.MakeCutList(outletface, [grains]) - outletshape = geompy.MakeScaleTransform(outletshape, oo, 1 / scale) - geompy.UnionList(outlet, geompy.SubShapeAll( - geompy.GetInPlace(shape, outletshape, True), geompy.ShapeType["FACE"])) - - symetry = [] - for (n, face) in enumerate(symetryface): - name = "symetry" + str(n) - symetry.append(geompy.CreateGroup(shape, geompy.ShapeType["FACE"], theName = name)) - symetryshape = geompy.MakeCutList(face, [grains]) - symetryshape = geompy.MakeScaleTransform(symetryshape, oo, 1 / scale) - geompy.UnionList(symetry[n], geompy.SubShapeAll( - geompy.GetInPlace(shape, symetryshape, True), geompy.ShapeType["FACE"])) + outlet = geompy.CreateGroup(shape, geompy.ShapeType["FACE"], theName = "outlet") + outletshape = geompy.MakeCutList(outletface, [grains]) + outletshape = geompy.MakeScaleTransform(outletshape, oo, 1 / scale) + geompy.UnionList(outlet, geompy.SubShapeAll( + geompy.GetInPlace(shape, outletshape, True), geompy.ShapeType["FACE"])) + + symetry = [] + for (n, face) in enumerate(symetryface): + name = "symetry" + str(n) + symetry.append(geompy.CreateGroup(shape, geompy.ShapeType["FACE"], theName = name)) + symetryshape = geompy.MakeCutList(face, [grains]) + symetryshape = geompy.MakeScaleTransform(symetryshape, oo, 1 / scale) + geompy.UnionList(symetry[n], geompy.SubShapeAll( + geompy.GetInPlace(shape, symetryshape, True), geompy.ShapeType["FACE"])) - groups = [] - groups.append(inlet) - groups.append(outlet) - groups.extend(symetry) + groups = [] + groups.append(inlet) + groups.append(outlet) + groups.extend(symetry) - if fillet: - strips = geompy.CreateGroup(shape, geompy.ShapeType["FACE"], theName = "strips") - shapeShell = geompy.ExtractShapes(shape, geompy.ShapeType["SHELL"], True) - stripsShape = geompy.MakeCutList(shapeShell[0], groups + [grainsOrigin]) - geompy.UnionList(strips, geompy.SubShapeAll( - geompy.GetInPlace(shape, stripsShape, True), geompy.ShapeType["FACE"])) - groups.append(strips) + if self.filletsEnabled: + strips = geompy.CreateGroup(shape, geompy.ShapeType["FACE"], theName = "strips") + shapeShell = geompy.ExtractShapes(shape, geompy.ShapeType["SHELL"], True) + stripsShape = geompy.MakeCutList(shapeShell[0], groups + [grainsOrigin]) + geompy.UnionList(strips, geompy.SubShapeAll( + geompy.GetInPlace(shape, stripsShape, True), geompy.ShapeType["FACE"])) + groups.append(strips) - wall = geompy.CutListOfGroups([sall], groups, theName = "wall") - groups.append(wall) + wall = geompy.CutListOfGroups([sall], groups, theName = "wall") + groups.append(wall) - return shape, groups + return shape, groups -def faceCentered(theta, fillet, direction): - if direction in [[1, 0, 0], [0, 0, 1]]: - return faceCenteredCubic(theta, fillet, direction) - - elif direction == [1, 1, 1]: - return faceCenteredHexagonalPrism(theta, fillet, direction) - - else: - raise Exception("This direction is not implemented") - diff --git a/anisotropy/genmesh.py b/anisotropy/genmesh.py index d3e4b9e..afb0e8c 100644 --- a/anisotropy/genmesh.py +++ b/anisotropy/genmesh.py @@ -72,12 +72,8 @@ def genmesh(root, name, direction, theta): # Shape ## geompy = getGeom() - structure = globals().get(p["name"]) - shape, groups = structure( - p["geometry"]["theta"], - p["geometry"]["fillets"], - p["geometry"]["direction"] - ) + structure = locals().get(p["name"]) + shape, groups = structure(**p["geometry"]) [length, surfaceArea, volume] = geompy.BasicProperties(shape, theTolerance = 1e-06) @@ -144,11 +140,17 @@ def genmesh(root, name, direction, theta): mesh.exportUNV(os.path.join(p["path"], "mesh.unv")) - stats = "" - for k, v in mesh.stats().items(): - stats += f"{ k }:\t\t{ v }\n" + meshStats = mesh.stats() + p["meshResult"] = dict( + mesh_id = p["mesh"]["id"], + surfaceArea = surfaceArea, + volume = volume, + **meshStats + ) + model.updateDB() - logger.info(f"mesh stats:\n{ stats[ :-1] }") + statstr = "\n".join(map(lambda k, v: f"{ k }:\t{ v }", meshStats)) + logger.info(f"mesh stats:\n{ statsstr }") salome.salome_close() diff --git a/anisotropy/models.py b/anisotropy/models.py index 047a3fb..5af8b98 100644 --- a/anisotropy/models.py +++ b/anisotropy/models.py @@ -83,6 +83,17 @@ class SubMesh(Mesh): class MeshResult(BaseModel): mesh_id = ForeignKeyField(Mesh, backref = "meshresults") + + surfaceArea = FloatField(null = True) + volume = FloatField(null = True) + + elements = FloatField(null = True) + edges = FloatField(null = True) + faces = FloatField(null = True) + volumes = FloatField(null = True) + tetrahedrons = FloatField(null = True) + prisms = FloatField(null = True) + pyramids = FloatField(null = True) calculationTime = TimeField(null = True) diff --git a/anisotropy/simple.py b/anisotropy/simple.py index e55be8c..8a89ca6 100644 --- a/anisotropy/simple.py +++ b/anisotropy/simple.py @@ -1,25 +1,21 @@ -#import salome -#salome.salome_init() - -#import GEOM -#from salome.geom import geomBuilder -#geompy = geomBuilder.New() - from math import pi, sqrt class simple(object): - def __init__(self, direction, theta, r0, L, radius, filletsEnabled, fillets): + def __init__(self, **kwargs): - self.direction = direction - self.theta = theta - self.r0 = r0 - self.L = L - self.radius = radius - self.filletsEnabled = filletsEnabled - self.fillets = fillets + self.direction = kwargs.get("direction", [1, 0, 0]) + self.theta = kwargs.get("theta", 0.01) + self.r0 = kwargs.get("r0", 1) + self.L = kwargs.get("L", 2 * self.r0) + self.radius = kwargs.get("radius", self.r0 / (1 - self.theta)) + self.filletsEnabled = kwargs.get("filletsEnabled", False) + self.fillets = kwargs.get("fillets", 0) def build(self): + from salomepl import getGeom + + geompy = getGeom() ### # Pore Cell @@ -164,7 +160,7 @@ class simple(object): # Groups ## shape = geompy.MakeCutList(poreCell, [grains]) - shape = geompy.MakeScaleTransform(shape, oo, 1 / scale, theName = "simpleCubic") + shape = geompy.MakeScaleTransform(shape, oo, 1 / scale, theName = "simple") sall = geompy.CreateGroup(shape, geompy.ShapeType["FACE"]) geompy.UnionIDs(sall,