De-la-mo

delamo package

Submodules

delamo.api module

Class hierarchy

  • DelamoModeler stores the various behind-the-scenes information required for model building
  • Part represents a CAD and finite element part
  • LayerPart is a particular kind of part
  • Assembly is a collection of parts and/or subassemblies
  • Layer is a particular type of Assembly, containing one or more LayerParts, representing a single layer
  • LaminaContact defines the contact/adhesion/cohesion between two lamina layeers
  • shell_solid_coupling stores the parameters needed for shell-to-solid coupling
  • CoordSys is an abstract class representing a coordinate system
  • SimpleCoordSys is an implementation of CoordSys that is a fixed Cartesian frame
class delamo.api.Assembly(**kwargs)

Bases: object

The assembly class represents an assembly of parts and assemblies.

classmethod FromAssemblies(name, *args)

Creates an assembly given the name parameter, and additional parameters representing the subassemblies to include in the assembly

classmethod FromParts(name, *args)

Creates an assembly given the name parameter, and additional parameters representing the parts to include in the assembly

classmethod FromPartsAndAssemblies(name, parts, assemblies)

Creates an assembly given the name parameter, and lists of parts and subassemblies to include in the assembly

MeshCohesive(*args, **kwargs)

This method iterates over the MeshCohesive methods of each part that is a child of this assembly (but not to subassemblies)

MeshSimple(*args, **kwargs)

This method iterates over the MeshSimple methods of each part that is a child of this assembly (but not to subassemblies)

assemblies = None

An ordered dictionary by name of Assembly objects or Assembly subclasses

name = None

The name of the assembly

partlist

This property gives a list of the parts in this assembly and all subassemblies

parts = None

An ordered dictionary by name of Part objects or Part subclasses

singlepart

Assuming this assembly contains only a single part, this property gives that part. Otherwise raises IndexError().

class delamo.api.CoordSys

Bases: object

Abstract class representing a coordinate system. Each CoordSys subclass should implement the method ApplyLayup(self,layerpart,layerdirection)

class delamo.api.DelamoModeler(**kwargs)

Bases: object

Represents the various behind-the-scenes information needed by the model builder

BodyNumDB = None

Actual body number database: dictionary indexed by body name of body number

BodyNumDB_build = None

ModelBuilder BodyNumber database proxy used in build phase

FEModel = None

Variable representing abq.mdb.models[‘FEModel’] which should be assigned in your initialization script

Finalize(script_to_generate, cad_file_path_from_script)

Write out the generated ABAQUS script and generated CAD file given the desired output filenames.

InitDict = None

A proxy object for a Python dictionary that is configured for execution with initinstrs

classmethod Initialize(globals, pointtolerance=1e-07, pointtolerancefactor=100.0, normaltolerance=0.01, tangenttolerance=0.01, GapWidth=0.3, STLMeshSize=3.0, Debug=False, license_key='')

Initialize the DelamoModeler, creating the codestores for initinstrs, assemblyinstrs, bcinstrs, meshinstrs, and runinstrs, and specifying the default tolerances

LaminateAssembly = None

Variable representing abq.mdb.models[‘FEModel’].rootAssembly which should be assigned in your initialization script

abaqus_assembly_functions(filename, globals)

Convenience routine for running a script in initinstrs Methods or attribute accesses for variables defined in this script will be referenced in accessinstrs (for the moment this is identical to abaqus_init_script()

abaqus_bc_functions(filename, globals)

Convenience routine for running a script in initinstrs Methods or attribute accesses for variables defined in this script will be referenced in bcinstrs

abaqus_fiber_functions(filename, globals)

Convenience routine for running a script in initinstrs Methods or attribute accesses for variables defined in this script will be referenced in fiberinstrs

abaqus_functions(globals)

Add the delamo builtin ABAQUS scripts abqfuncs_assembly.py, abqfuncs_bc.py, abqfuncs_mesh.py, and abqfuncs_run.py to the DelamoModeler codestores

abaqus_init_script(filename, globals)

Convenience routine for running a script in initinstrs Methods or attribute accesses for variables defined in this script will be referenced in assemblyinstrs

abaqus_mesh_functions(filename, globals)

Convenience routine for running a script in initinstrs Methods or attribute accesses for variables defined in this script will be referenced in meshinstrs

abaqus_run_functions(filename, globals)

Convenience routine for running a script in initinstrs Methods or attribute accesses for variables defined in this script will be referenced in runinstrs

abq = None

Proxy object represents the “abq” module, for execution with initinstrs

abq_assembly = None

Proxy object represents the “abq” module, for execution with assemblyinstrs

abqpointtolerance = None

Default point tolerance for finite element face and edge identification

add_abaqus_code(instrs, func_code, globals)

Add the specified ABAQUS code func_code to the specified instrs codestore. Proxies for the left sides of any assignment statements found in func_code will be added to the dictionary globals

assembly = None

Proxy object represents the ABAQUS “assembly” module, for execution with assemblyinstrs

assemblyinstrs = None

Codestore of ABAQUS assembly instructions

autofiber = None

Abaqus reference to the auto fiber class

bcinstrs = None

Codestore of ABAQUS boundary condition instructions

connector = None

Proxy object represents the ABAQUS “connector” module, for execution with bcinstrs

exec_abaqus_script(execinstrs, accessinstrs, filename, globals)

Add the ABAQUS code from the specified file to the specified execinstrs codestore. Proxies for the left sides of any assignment statements found in func_code will be added to the dictionary globals, for access in the accessinstrs codestore.

fiberinstrs = None

Codestore of ABAQUS fiber orientation instructions

get_unique()

Get a unique number. The numbers from this method of a particular DelamoModeler object will never repeat.

globals = None

Global variable dictionary providng context for main script and imported functions/methods, including functions from abqfuncs_*.py

initinstrs = None

Codestore of ABAQUS initialization instructions

mesh = None

Proxy object represents the ABAQUS “mesh” module, for execution with meshinstrs

meshinstrs = None

Codestore of ABAQUS meshing instructions

modelbuilder = None

delamo.OCCModelBuilder.OCCModelBuilder() object.

normaltolerance = None

Default normal tolerance for finite element face identification

regionToolset = None

Proxy object represents the ABAQUS regionToolset module, for execution with assemblyinstrs

runinstrs = None

Codestore of ABAQUS running instructions

section = None

Proxy object represents the ABAQUS “section” module, for execution with initinstrs

stepfile = None

Variable wrapped for ABAQUS access representing name of generated CAD file to load

tangenttolerance = None

Default tangent tolerance for finite element edge identification

to_be_saved = None

List of MBBody objects to be saved into the output CAD file

to_be_saved_dependencies = None

MBBody objects in to_be_saved may be dependent on the existence of their parent objects. This is a list, a place to store such references so their parents don’t get destroyed

uniquenumber = None

value for next unique number from get_unique()

class delamo.api.LaminaContact(**kwargs)

Bases: object

LaminaContact represents the Contact/Adhesion/Cohesion between two layers. It includes storage of the various tie/contact/cohesive objects tying the two layers together

DM = None

Reference to the DelamoModeler object

DefineContinuity(FEModel, bottom_laminapart, top_laminapart, point_normal, pointtolerance, normaltolerance, master=None)

Define a TIE (continuity) boundary condition between a face of a part of the bottom lamina and a corresponding face of a part of the top lamina

DefineDelamination(FEModel, ContactInteraction, bottom_laminapart, top_laminapart, point_normal, pointtolerance, normaltolerance, master=None)

Define a contact boundary condition between a face of a part of the bottom lamina and a corresponding face of a part of the top lamina

DefineLamination(FEModel, CohesiveInteraction, bottom_laminapart, top_laminapart, point_normal, pointtolerance, normaltolerance, master=None)

Define a cohesive boundary condition between a face of a part of the bottom lamina and a corresponding face of a part of the top lamina

bottomlamina = None

Bottom layer

cohesives = None

list of contact objects attaching these two layers

contacts = None

list of contact objects between these two layers

ties = None

list of Tie objects attaching these two layers

toplamina = None

Top layer

class delamo.api.Layer(**kwargs)

Bases: delamo.api.Assembly

Represents a layer of a composite material

CreateFiberObject(DM, point, fibervec, normal, mp, fiberint=1.0, angle_error=0.01, final_plotting=False)

Utilize Autofiber orientation package to calculate optimal fiber orientations at each mesh element centroid.

classmethod CreateFromMold(DM, mold, direction, thickness, name, Section, layup, meshsize=0.5, coordsys=None)

Create a layer atop the specified mold. * direction: “OFFSET” or “ORIG” * thickness: Thickness of layer (offsetting operation) * name: Unique name for layer * Section: ABAQUS section fo the layer * layup: Ply orientation in degrees * coordsys: Reference coordinate system for layup

classmethod CreateFromParams(DM, create_params, name, LayerSection, layupdirection, meshsize=0.5, split=None, coordsys=None)

Create a layer given creation parameters to be passed to the geometry kernel, a name, ABAQUS Section, layup direction, etc.

Finalize(DM)

Build Python and Finite Element structure for this layer based on geometry kernel structure. It is not permissible to break a layer into pieces (multiple parts or bodies) after this finalize call. It IS permissible to do surface imprints that operate in-place on the existing layerbodies.

GetPartInstanceFaceRegionFromPoint(Point, PointTolerance)
LayerSection = None

ABAQUS Section used by this layer

LayupFiberObject(DM, layupdirection, final_plotting=False)

Use Autofiber object to determine layup orientation based on a fiber angle.

Split(splitpath_filename, PointTolerance)
coordsys = None

Coordinate system used by this layer

gk_layer = None

Underlying geometry kernel object

layupdirection = None

Orientation of this layer, in degrees

class delamo.api.LayerPart(**kwargs)

Bases: delamo.api.Part

A LayerPart is a Part that will be used in a Layer

ApplyLayup(coordsys, layupdirection, fiberorientation)

Assign self.fe_datum_csys and self.fe_materialorientation

class delamo.api.Part(**kwargs)

Bases: object

The Part class represents a CAD part and a parallel finite element (ABAQUS) part.

AssignSection(Section, offset=0.0, offsetType=<delamo.genericobjectwrapper.wrappedobject object>)

Assign ABAQUS finite element Section (which must be an Abaqus section object with a material applied to it) to this part. Calls SectionAssigment() method on the part

CreateInstance(dependent=<delamo.genericobjectwrapper.wrappedobject object>)

Create the finite element instance of this part

DM = None

DelamoModeler object

classmethod FromGK3D(DM, gk_layerbody, shell=False, no_FE_instance=False, omit_from_FE=False)

Define a solid or shell part from a 3D geometry kernel object. Ensures that this geometry kernel object will be saved to the generated CAD file and programs ABAQUS to load it in and create an instance.

classmethod FromGeometryFile(DM, filename, name, bodynum)

Define a finite element solid or shell part from a file on disk. This will not have a geometry kernel representation

GetInstanceEdge(edgepoints, pointtolerance)

Get an instance edge based on a point and tolerance.

GetInstanceEdgeRegion(edgepoints, pointtolerance)

Get an instance edge region based on a point and tolerance. Note that this creates a “Set-like” region, not the “Surface-like” regions used in the bonding functions(below) in api.py. See Abaqus Scripting Reference Guide section 45.3 for the distinction

GetInstanceEdgeRegionSurface(edgepoints, pointtolerance)

Get an instance edge region based on a point and tolerance. Note that this creates a “Surface-like” region with side1Edges= used in the bonding functions(below) in api.py, not the “Set-like” regions used in some other contexts. See Abaqus Scripting Reference Guide section 45.3 for the distinction.

GetInstanceEdgeRegion_point_tangent(edgepointtangent, pointtolerance, tangenttolerance)

Get an instance edge region based on a point and tangent and tolerances. Note that this creates a “Set-like” region, not the “Surface-like” regions used in the bonding functions(below) in api.py. See Abaqus Scripting Reference Guide section 45.3 for the distinction

GetInstanceEdge_point_tangent(edgepointtangents, pointtolerance, tangenttolerance)

Get an instance edge based on a point and tangent and tolerances.

GetInstanceFace(surface_points, pointtolerance)

Get an instance face based on a point and tolerance.

GetInstanceFaceRegion(surface_points, pointtolerance)

Get an instance face region based on a point and tolerance. Note that this creates a “Set-like” region, not the “Surface-like” regions used in the bonding functions(below) in api.py. See Abaqus Scripting Reference Guide section 45.3 for the distinction

GetInstanceFaceRegionSurface(surface_points, pointtolerance)

Get an instance face region based on a point and tolerance. Note that this creates a “Surface-like” region used in the bonding functions(below) in api.py, not the “Set-like” regions used in some other contexts. See Abaqus Scripting Reference Guide section 45.3 for the distinction.

GetInstanceFaceRegion_point_normal(surface_point_and_normal, pointtolerance, normaltolerance)

Get an instance face region based on a point and normal and tolerances. Note that this creates a “Set-like” region, not the “Surface-like” regions used in the bonding functions(below) in api.py. See Abaqus Scripting Reference Guide section 45.3 for the distinction

GetInstanceFace_point_normal(surface_point_and_normal, pointtolerance, normaltolerance)

Get an instance face based on a point and normal and tolerances.

GetMultipleInstanceEdges(edgepoints, pointtolerance)

Get multiple instance edges based on points and tolerance.

GetMultipleInstanceEdgesRegion(edgepoints, pointtolerance)

Get instance edge regions based on a point and tolerance. Note that this creates a “Set-like” region, not the “Surface-like” regions used in the bonding functions(below) in api.py. See Abaqus Scripting Reference Guide section 45.3 for the distinction

GetMultipleInstanceFaces(surface_points, pointtolerance)

Get instance faces based on points and a tolerance.

GetMultipleInstanceFacesRegion(surface_points, pointtolerance)

Get an instance face region based on multile faces selected by surface points and tolerance. Note that this creates a “Set-like” region, not the “Surface-like” regions used in the bonding functions(below) in api.py. See Abaqus Scripting Reference Guide section 45.3 for the distinction

GetMultiplePartEdges(edgepoints, pointtolerance)

Method to find multiple edges of the part given a point and tolerance

GetMultiplePartFaces(surface_points, pointtolerance)

Method to find multiple faces of the part given a point and tolerance

GetPartEdge(edgepoints, pointtolerance)

Method to find an edge of the part given a point and tolerance

GetPartEdge_point_tangent(edgepointtangent, pointtolerance, tangenttolerance)

Method to find an edge of the part given a point and tangent and tolerances

GetPartFace(surface_points, pointtolerance)

Method to find a face of the part given a point and tolerance

GetPartFace_point_normal(surface_point_and_normal, pointtolerance, normaltolerance)

Method to find a face of the part given a point and normal and tolerances

MeshCohesive(meshsize, ElemShape=None, Algorithm=None, ElemLibrary=None, refined_edges=[], pointtolerance=None, refinedmeshsize=None, DeviationFactor=None, MinSizeFactor=None, SweepSense=None)

Meshing routine for cohesive layers: meshsize is nominal meshing size ElemShape should be abqC.HEX or abqC.HEX_DOMINATED (default) Algorithm can be abqC.ADVANCING_FRONT (default) or abqC.MEDIAL_AXIS ElemLibrary can be abqC.STANDARD (default) or abqC.EXPLICIT refined edges is a list of tuples: (endpoint1, interiorpoint, endpoint2) representing 3 points on each edge to be specially refined pointtolerance is tolerance size for finding refined edges refinedmeshsize is size of refined mesh regions DeviationFactor and MinSizeFactor are scaling factors for seeding (default 0.1) SweepSense should be abqC.FORWARD (default) or abqC.REVERSE

MeshSimple(ElemTypes, meshsize, ElemShape=None, ElemTechnique=None, refined_edges=[], pointtolerance=None, tangenttolerance=None, refinedmeshsize=None, DeviationFactor=None, MinSizeFactor=None)

Perform meshing of this Part. meshsize is nominal meshing size ElemShape should be abqC.HEX, abqC.HEX_DOMINATED (default), abqC.TET, etc. for solids or abqC.QUAD or abqC.QUAD_DOMINATED (default), or abqC.TRI for shells ElemTechnique can be abqC.SYSTEM_ASSIGN (default), abqC.FREE, or abqC.STRUCTURED refined edges is a list of tuples: ((point1, tangent1),(point2,tangent2)) on each edge to be specially refined pointtolerance is tolerance size for finding refined edges refinedmeshsize is size of refined mesh regions DeviationFactor and MinSizeFactor are scaling factors for seeding (default 0.1)

SeedPartEdgesByFaces(surface_points_and_normals, pointtolerance, normaltolerance, meshsize)

Seed the edges around a set of faces with a particular meshing size. This is used for localized mesh refinement. The actual implementation is the ABAQUS code in abqfuncs_mesh.py

fe_datum_csys = None

Abaqus datum CSYS

fe_inst = None

Wrapper for Abaqus part instance, assembly phase

fe_materialorientation = None

Abaqus material orientation

fe_part = None

Wrapper for ABAQUS part object, assembly phase

fe_part_meshing = None

Wrapper for Abaqus part object, meshing phase

gk_layerbody = None

Geometry kernel representation of this part

name = None

Name of this part

shell = None

True/False: Is a shell (as opposed to a solid) part

class delamo.api.SimpleCoordSys(fibervec, crossfibervec)

Bases: delamo.api.CoordSys

Concrete implementation of a CoordSys representing a fixed Cartesian coordinate frame

ApplyLayup(layerpart, layupdirection)

Set the material orientation of the specified layerpart to layupdirection (in degrees) relative to this coordinate frame

crossfibervec = None

Unit vector along the fibers of a 90 degree ply

fibervec = None

Unit vector along the fibers of a 0 degree ply

outofplanevec = None

Out-of-plane unit vector, i.e. fibervec cross crossfibervec

delamo.api.bond_layers(DM, layer1, layer2, defaultBC='TIE', delamBC='CONTACT', delamRingBC='NONE', CohesiveInteraction=None, ContactInteraction=None, delaminationlist=None, master_layer=None, cohesive_layer=None, delamo_sourceline=None, delamo_phase=None, delamo_basename=None)

Bond two layers together. Parameters: * DM: DelamoModeler object * layer1: First layer * layer2: Second layer * defaultBC: The boundary condition for the bonded zone: Generally “TIE”, “COHESIVE”, or “COHESIVE_LAYER” * delamBC: The boundary condition for the bulk of the delaminated region(s). Generally “CONTACT” * delamRingBC: The boundary condition for the outer zone of the delaminated region: Generally “NONE” * CohesiveInteraction: The ABAQUS interaction property for any cohesive portions of the bond * ContactInteraction: The ABAQUS interaction property for any contact portions of the bond * delaminationlist: A list of files with delamination outlines (loops of 3D coordinates) * master_layer: A preference for which layer should be the master layer in ABAQUS * cohesive_layer: If defaultBC==”COHESIVE_LAYER”, this layer will be used to achieve the bond. The cohesive_layer should have finite thickness but should NOT be finalized (this routine will finalize the cohesive_layer). You will also need to mesh the cohesive_layer afterward using MeshCohesive() * delamo_sourceline: An index of the line number of this bond_layers() call, used to distinguish between mulitiple bonding steps. * delamo_phase: Phase of this multi-step defect insertion process. * delamo_basename: Base directory name for creation of layer boundary (STL) meshes.

delamo.api.shell_and_cutout_from_shelltool(DM, shelltool_filename)

Load in a .SAT file with two CAD objects. Both are presumed to be sheet bodies (shells) They should both be colored. Exactly one of them should have a red component greater than 0.75. The one without a red component greater than 0.75 is the “shell”; the one with a red component greater than 0.75 is the “tool”, which is assumed to intersect with the shell, yielding a closed curve. This function performs the intersection of shell and tool, and makes a shell with a hole, plus a cut-out that fits in the hole. This function returns a Part reprepresenting the shell with hole, a second Part representing the cut-out, and a list of (point,tangent) tuples for identifying the edges around the hole (useful for mesh refinement)

class delamo.api.shell_solid_coupling(**kwargs)

Bases: object

This class stores the various parameters required to implement shell-to-solid coupling

Finalize(DM, influenceDistance=None, positionTolerance=None)

Create the ABAQUS ShellSolidCoupling object.

bond_layer(DM, layer, layeroffset)

Bond a layer to the shell

cutout = None

Cutout representation, the second part returned from shell_and_cutout_from_shelltool()

classmethod from_shell_and_cutout(DM, shell, cutout)

Create shell_solid_coupling object from a shell and cutout, for example as returned by shell_and_cutout_from_shelltool()

shell = None

Shell representation with hole, the first Part returned from shell_and_cutout_from_shelltool()

shell_edge_normal_list = None

List of normals on the edges of the shell cutout (in the plane of the shell, will be normal to the edge of the layer).

shell_edge_point_list = None

List of points on the edges of the shell cutout

shell_edge_tangent_list = None

List of tangents on the edges of the shell cutout

surfacefaces = None

List of lists of ABAQUS side1faces type face sets, from all shell/cutout edges

delamo.simple module

delamo.simple.laminate(initscript, controlpoints, layup, fixedpoint=None, forcesurfacepoint=None, delamination=None, delam_index=0)

Module contents