输出平均应力应变

EDC21

我一个模型是15*10*10的长方体，有1500个单元。有没有大神能写一个代码来遍历这个长方体内所有的n*n*n的小正方体（从1*1*1到10*10*10）比如1*1*1的有1500,2*2*2的有1134个，3*3*3的有832个,然后输出每个正方体区域的平均应力应变（Stress average = Sum (Stress (element)*volume (element))/Sum (volume (element))）。
例如一个2*2*2的正方体， 拿由单元编号1，2，11，12，101，102，111，112组成的2*2*2的来说，我要得到这八个单元的应力，分别乘以他们的单元体积，然后除以他们的总体积.
有一份参考代码
import sys, getopt, os, string
import math
from odbAccess import *
from abaqusConstants import *

odbPath = "D:\ho\Desktop\Abaqus Script\M1.odb"  
odb = session.openOdb(name=odbPath,readOnly=FALSE)      
# retrieve steps from the odb
grout_instance = odb.rootAssembly.instances['COULIS_SMOOTH-1']
numElem = len(grout_instance.elements)
keys = odb.steps.keys()
for stp in keys:
        step = odb.steps[stp]
        # retrieve frames from the odb
        frameRepository = step.frames
        numFrames = len(frameRepository)
        for fr in range(0,numFrames):
                frame=step.frames[fr]
                print 'Id = %d, Time = %f\n'%(frame.frameId,frame.frameValue)
                # get fieldOutputs object
                fo = frame.fieldOutputs
                S  = fo['S']
                EVOL = fo['EVOL']
                S_grout = S.getSubset(region=grout_instance,\
                                                position=INTEGRATION_POINT,
                                                elementType='CAX3')
                EVOL_grout = EVOL.getSubset(region=grout_instance)
                myStress_data = {} #Declaration of my New stress Data
                # Loops over elements to get stpress S11 and Volume
                for i in range(0,len(S_grout.values)):
                        evol= EVOL_grout.values.data
                        a = (S_grout.values.data[0])*evol
                        elemId=S_grout.values.elementLabel
                        myStress_data.setdefault(elemId,[]).append(a)
                # Append
                elementLabels = []
                elementData = []
                for key in sorted(myStress_data.iterkeys()):
                        elementLabels.append(key)
                        elementData.append(myStress_data[key])
                myS11 = frame.FieldOutput(name='myS11',
                                                                           description='my Principle stpress S11',
                                                                           type=SCALAR)
                myS11.addData(position=INTEGRATION_POINT,
                                                   instance=grout_instance,
                                                   labels=elementLabels,
                                                   data=elementData)
odb.save()
odb.close()

maschinen666

建模设置各自对应的组， 然后输出每单个组的应力结果就好了, 在dat file里面可已找到每一个对应单元的应力结果啊.