Abaqus Python Script

A friend asked me to write a python script to control the FEM software Abaqus. He wanted to simulate a compound material consisting of a material block with small spheres of a diffrent material inserted. It was the first time I ever saw Abaqus, but got the hang of it quite fast.

The resulting script produces a block of material, generates n spheres wich are distributed randomly and not overlapping inside the block. The radius of the spheres is controled by their number. It is choosen in a way that the n spheres make up 20% of the volume of the block.
Sphere Positions and radii are written out into a file. Im looking forward to his simulation results.

from part import *
from material import *
from section import *
from assembly import *
from step import *
from interaction import *
from load import *
from mesh import *
from job import *
from sketch import *
from visualization import *
from connectorBehavior import *
import random
import math

#Euclidean Distance Measure
def euclidean_py(x, y):
    # lightly modified from implementation by Thomas Sicheritz-Ponten.
    # This works faster than the Numeric implementation on shorter
    # vectors.
    if len(x) != len(y):
        raise ValueError, "vectors must be same length"
    sum = 0
    for i in range(len(x)):
        sum += (x[i]-y[i])**2
    return math.sqrt(sum)

#Create Cube
mdb.models['Model-1'].ConstrainedSketch(name='__profile__', sheetSize=200.0)
mdb.models['Model-1'].sketches['__profile__'].rectangle(point1=(0.0, 0.0), point2=(1.0, 1.0))
mdb.models['Model-1'].Part(dimensionality=THREE_D, name='Part-1', type=DEFORMABLE_BODY)
mdb.models['Model-1'].parts['Part-1'].BaseSolidExtrude(depth=1.0, sketch=mdb.models['Model-1'].sketches['__profile__'])
del mdb.models['Model-1'].sketches['__profile__']

#Instantiate cube
mdb.models['Model-1'].rootAssembly.Instance(dependent=ON, name='tmpCube', part=mdb.models['Model-1'].parts['Part-1'])

#File Output
text_file = open("c:\Sphere_Positions.txt", "w")

#Number of Spheres
n = 10
text_file.write("Number of Spheres: " + str(n) + "\n")

#Radius of one sphere
r = (3*0.20)/(4*math.pi*n)
r = r**(1.0/3.0)
print "Sphere Radius is " + str(r)
text_file.write("Sphere Radius is " + str(r) + "\n")

#Create Sphere
mdb.models['Model-1'].ConstrainedSketch(name='__profile__', sheetSize=200.0)
mdb.models['Model-1'].sketches['__profile__'].ConstructionLine(point1=(0.0, -100.0), point2=(0.0, 100.0))
mdb.models['Model-1'].sketches['__profile__'].FixedConstraint(entity=mdb.models['Model-1'].sketches['__profile__'].geometry[2])
mdb.models['Model-1'].sketches['__profile__'].ArcByCenterEnds(center=(0.0, 0.0), direction=CLOCKWISE, point1=(0.0, -r), point2=(0.0, r))
mdb.models['Model-1'].sketches['__profile__'].Line(point1=(0.0, -r), point2=(0.0, r))
mdb.models['Model-1'].sketches['__profile__'].VerticalConstraint(entity=mdb.models['Model-1'].sketches['__profile__'].geometry[4])
mdb.models['Model-1'].Part(dimensionality=THREE_D, name='Part-2', type=DEFORMABLE_BODY)
mdb.models['Model-1'].parts['Part-2'].BaseSolidRevolve(angle=360.0, flipRevolveDirection=OFF, sketch=mdb.models['Model-1'].sketches['__profile__'])
del mdb.models['Model-1'].sketches['__profile__']

#List of Spheres
sphereList = []
sphereInstancesList = []

#Create n instances of the sphere
for i in range(1, n+1):
	InstanceName = 'Sphere_' + str(i)
	print InstanceName
	text_file.write(InstanceName)

	#Maximum tries to distribute sphere
	maxTries = 1000

	while len(sphereList) < i:
		maxTries -= 1
		if maxTries < 1:
			print "Maximum Distribution tries exceded. Error! Restart the Script!"
			break;

		#Make sure Spheres dont cut cube sides
		vecPosition = [r+(random.random()*(1.0-r-r)),r+(random.random()*(1.0-r-r)),r+(random.random()*(1.0-r-r))]
		for pos in sphereList:
			if euclidean_py(pos, vecPosition) < 2*r:
				break
		else:
			sphereList.append(vecPosition)
			print vecPosition
			text_file.write("\t" + str(vecPosition) + "\n")
			#Instantiate Sphere
			mdb.models['Model-1'].rootAssembly.Instance(dependent=ON, name=InstanceName, part=mdb.models['Model-1'].parts['Part-2'])
			#Translate Instance of Sphere
			mdb.models['Model-1'].rootAssembly.translate(instanceList=(InstanceName, ), vector=vecPosition)
			sphereInstancesList.append(mdb.models['Model-1'].rootAssembly.instances[InstanceName])

print "Sphere Radius is " + str(r)

#Cut all spheres with cube
mdb.models['Model-1'].rootAssembly.PartFromBooleanCut(cuttingInstances=(
sphereInstancesList ),
instanceToBeCut=mdb.models['Model-1'].rootAssembly.instances['tmpCube'],
name='Part-3')

#Intantiate Cut Object
mdb.models['Model-1'].rootAssembly.Instance(dependent=ON, name='Cube',
part=mdb.models['Model-1'].parts['Part-3'])

#Delete temporary cube
del mdb.models['Model-1'].rootAssembly.instances['tmpCube']

text_file.close()

Filed under: Abaqus, Python

Gray Hat Python (error 50 on 64 bit)

The great Book “Gray Hat Python – Python Programming for Hackers and Reverse Engineers” by Justin Seitz includes some great sample code to build your own debugger using only Python, ctypes and the win32 api. The exmple source code is freely available on the book website at nostarch press: http://nostarch.com/ghpython.htm

Unfortunatly the Windows Debugging Tools, namly the API functions don’t work for 64 bit versions of Windows. So be prepared if  the DebugActiveProcess call does return Null.

To check if the failed function call is caused by API problems insert the following just behind the “else” which indicates the failed attachment to the debug process:

print "[*] Unable to attach to the process. %s" % FormatError(kernel32.GetLastError())

The resulting error has the error code 50 (ERROR_NOT_SUPPORTED). The text of the message varies depending your OS language. In german it would be “Die Anforderung wird nicht unterstützt.”, in english “The request is not supported.”.

A solutions to get the examples running would be to use the Windows XP Mode bundled with Win 7 or any other virtual 32 bit machine.

Filed under: Python