The script can be used to apply periodic boundary conditions in 2D Recatangle
Representative Volume Element (RVE) of any size. One has to provide certain inputs,
1. model name of the RVE and 2. height and 3. width of RVE,
In order to use this script to your FEA model, you need to change the modelkey variable name to
a name of your model.
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 | # Applying Periodic Boundary Conditions in Abaqus/CAE # Created by: Paramveer Sharma # MS, Aerospace Engineering # Indian Institute of Technology Madras # PBC Constraints is generalised for any size of rve # You can change the model name : Model-1 to whatever the name of your model and same apply for Instance, ######################################### modelkey ='Model-Name' instancekey = mdb.models[modelkey].rootAssembly.instances.keys()[0] F=mdb.models[modelkey].rootAssembly.instances[instancekey].nodes cornsetlist,setlabel=['Bottemleft','Bottemright','Topright','Topleft'],['setlabelH','consH','setlabelV','consV','corn','U','V'] a,noj = getInput('length of RVE (x - coordinates)'),0 b = getInput('height of RVE (y - coordinates)') a= float(a) b = float(b) for i in range(len(F)): coor=F[i].coordinates column = [0,a,b] for k in column: if coor[0]==k: for m in column: if coor[1]==m: #print(str(F[i].coordinates)) cornseti=setlabel[4]+str(F[i].label) setsi=mdb.models[modelkey].rootAssembly.Set(cornseti, F[i:i+1]) if coor[0]==0.0 and coor[1]==0.0: cornsetlist[0]=cornseti elif coor[0]==a and coor[1]==0.0: cornsetlist[1]=cornseti elif coor[0]==a and coor[1]==b: cornsetlist[2]=cornseti elif coor[0]==0.0 and coor[1]==b: cornsetlist[3]=cornseti if noj==0: mdb.models[modelkey].Equation('corneq1',([1,cornsetlist[2],1],[-1,cornsetlist[1],1],[-1,cornsetlist[3],1],[1,cornsetlist[0],1])) mdb.models[modelkey].Equation('corneq2',([1,cornsetlist[3],2],[-1,cornsetlist[0],2],[-1,cornsetlist[2],2],[1,cornsetlist[1],2])) if noj==1: mdb.models[modelkey].Equation('corneq1',([1,cornsetlist[1],1],[-1,cornsetlist[2],1],[-1,cornsetlist[3],1],[1,cornsetlist[0],1])) mdb.models[modelkey].Equation('corneq2',([1,cornsetlist[2],2],[-1,cornsetlist[1],2],[-1,cornsetlist[3],2],[1,cornsetlist[0],2])) if noj==2: mdb.models[modelkey].Equation('corneq1',([1,cornsetlist[2],1],[-1,cornsetlist[1],1],[-1,cornsetlist[3],1],[1,cornsetlist[0],1])) mdb.models[modelkey].Equation('corneq2',([1,cornsetlist[2],2],[-1,cornsetlist[1],2],[-1,cornsetlist[3],2],[1,cornsetlist[0],2])) column = [0,a,b] l1,l2,l3,l4=[],[],[],[] ###l1,l2,l3,l4 is list of node index position which situated at Bottem, Top, Right and Left respectively. for i in range(len(F)): coor=F[i].coordinates if coor[1]==0: if coor[0]!=0 and coor[0]!=a: l2.append(i); if coor[1]==b: if coor[0]!=0 and coor[0]!=a: l1.append(i); if coor[0]==0: if coor[1]!=0 and coor[1]!=b: l4.append(i); if coor[0]==a: if coor[1]!=0 and coor[1]!=b: l3.append(i); for i in range(len(l1)): coor1=F[l1[i]].coordinates for j in range(len(l2)): coor2=F[l2[j]].coordinates if coor1[0] == coor2[0]: setlabeli=setlabel[0]+str(F[l1[i]].label) setlabelj=setlabel[0]+str(F[l2[j]].label) cequation1=setlabel[1]+str(F[l1[i]].label)+setlabel[5]+str(F[l2[j]].label) cequation2=setlabel[1]+str(F[l1[i]].label)+setlabel[6]+str(F[l2[j]].label) setsi=mdb.models[modelkey].rootAssembly.Set(setlabeli, F[l1[i]:l1[i]+1]) setsj=mdb.models[modelkey].rootAssembly.Set(setlabelj, F[l2[j]:l2[j]+1]) mdb.models[modelkey].Equation(cequation1,([1,setlabeli,1],[-1,setlabelj,1],[-1,cornsetlist[3],1],[1,cornsetlist[0],1])) mdb.models[modelkey].Equation(cequation2,([1,setlabeli,2],[-1,setlabelj,2],[-1,cornsetlist[3],2],[1,cornsetlist[0],2])) for i in range(len(l3)): coor1=F[l3[i]].coordinates for j in range(len(l4)): coor2=F[l4[j]].coordinates if coor1[1] == coor2[1]: setlabeli=setlabel[2]+str(F[l3[i]].label) setlabelj=setlabel[2]+str(F[l4[j]].label) cequation1=setlabel[3]+str(F[l3[i]].label)+setlabel[5]+str(F[l4[j]].label) cequation2=setlabel[3]+str(F[l3[i]].label)+setlabel[6]+str(F[l4[j]].label) setsi=mdb.models[modelkey].rootAssembly.Set(setlabeli, F[l3[i]:l3[i]+1]) setsj=mdb.models[modelkey].rootAssembly.Set(setlabelj, F[l4[j]:l4[j]+1]) mdb.models[modelkey].Equation(cequation1,([1,setlabeli,1],[-1,setlabelj,1],[-1,cornsetlist[1],1],[1,cornsetlist[0],1])) mdb.models[modelkey].Equation(cequation2,([1,setlabeli,2],[-1,setlabelj,2],[-1,cornsetlist[1],2],[1,cornsetlist[0],2])) #################################################################### #PBC Constraints end ################################################################### |
This script can be used to extract the nodel point information such as node label, X-Y coordinates of nodes
jacobian of isoparametric transformation of elements deformation. All the information extracted can be saved into the text file generated by
this script
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 | #_________________________________________________________ #_____ Script for finding the nodal connectivity in abaqus #_________________________________________________________ # Paramveer Sharma # MS, Aerospace Engineering # Indian Institute of Technology Madras # http://pyparam.github.io # Email: paramveersharma9@gmail.com import part import numpy as np for ii in range(4): liststring = ['Glabel', 'XCoor','YCoor', 'Jacobian'] file = open('/home/paramveer/'+liststring[ii]+'.txt','w+') # path and filename you want to generate elemArr = mdb.models['Model-1-2-1'].parts['Part-1'].sets['setsection2'].elements # pick the part you want to analyse len(elemArr) print len(elemArr) counter = 0 for e in elemArr: ktuple = e.getNodes() length = len(ktuple) elabel = e.label counter = counter +1 if length ==4: klabel1 = ktuple[0].label kcoor1 = ktuple[0].coordinates klabel2 = ktuple[1].label kcoor2 = ktuple[1].coordinates klabel3 = ktuple[2].label kcoor3 = ktuple[2].coordinates klabel4 = ktuple[3].label kcoor4 = ktuple[3].coordinates s, t= 0 , 0 X = np.array([[kcoor1[0], kcoor2[0], kcoor3[0], kcoor4[0]]]) Y = np.transpose([[kcoor1[1],kcoor2[1],kcoor3[1],kcoor4[1]]]) guass_x = (kcoor1[0] + kcoor2[0] +kcoor3[0] + kcoor4[0])/4.0 guass_y = (kcoor1[1] + kcoor2[1] +kcoor3[1] + kcoor4[1])/4.0 r1 , r2 = [0, 1-t, t-s, s-1], [t-1, 0 , s+1, -s-t] r3 , r4 = [s-t, -s-1, 0, t+1], [1-s, s+t, -t-1, 0] isomatrix = np.array([r1,r2,r3,r4]) temp = X.dot(isomatrix) jacobian = (1.0/8.0)*(temp.dot(Y)) area = (1.0/2.0)*(kcoor1[0]*kcoor2[1]+kcoor2[0]*kcoor3[1]+kcoor3[0]*kcoor4[1]+kcoor4[0]*kcoor1[1]-kcoor2[0]*kcoor1[1]-kcoor3[0]*kcoor2[1]-kcoor4[0]*kcoor3[1]-kcoor1[0]*kcoor4[1]) if area<0.0: area = (-1.0)*area listout = [elabel, guass_x,guass_y,jacobian] file.write('%g\n' % (listout[ii])) else: klabel1 = ktuple[0].label kcoor1 = ktuple[0].coordinates klabel2 = ktuple[1].label kcoor2 = ktuple[1].coordinates klabel3 = ktuple[2].label kcoor3 = ktuple[2].coordinates zeta1, zeta2, zeta3= (1.0/3.0) , (1.0/3.0), (1.0/3.0) X = np.array([[kcoor1[0], kcoor2[0], kcoor3[0]]]) Y = np.transpose([[kcoor1[1],kcoor2[1],kcoor3[1]]]) guass_x = kcoor1[0]*zeta1 +kcoor2[0]*zeta2 +kcoor3[0]*zeta3 guass_y = kcoor1[1]*zeta1 +kcoor2[1]*zeta2 +kcoor3[1]*zeta3 jacobian = (1.0/2.0)*((kcoor2[0]*kcoor3[1]-kcoor2[1]*kcoor3[0])+(kcoor3[0]*kcoor1[1]-kcoor1[0]*kcoor3[1]) +(kcoor1[0]*kcoor2[1]-kcoor2[0]*kcoor1[1])) area = jacobian listout = [elabel, guass_x,guass_y,jacobian] file.write('%g\n' % (listout[ii])) file.close() |