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