update case and odb2vtk

.gitignore

@@ -0,0 +1,12 @@
+*.pyc
+*.rpy*
+*.msg
+*.com
+*.dat
+*.ipm
+*.jnl
+*.lck
+*.log
+*.prt
+*.sim
+*.sta

source/odb2vtk-master/README.md

@@ -0,0 +1,5 @@
+# odb2vtk
+Python script converts Abaqus ODB files to VTK format for Paraview visualisation. 
+This Python script is released with the paper:
+Qingbin Liu, Jiang Li, Jie Liu, 2017, ParaView visualization of Abaqus output on the mechanical deformation of complex microstructures, Computers and Geosciences, 99: 135-144.
+Please cite it if the script is used in your reports/papers.

source/odb2vtk-master/odb2vtk.txt

@@ -0,0 +1,12 @@
+----------input and output path----------
+odb_path = 'C:\Temp'
+odb_name = 'CP10_L6_DP1'
+vtk_path = 'C:\Temp\CP10_L6_DP1'
+--------------type of mesh--------------
+mesh_type = '12'
+-------------number of piece-------------
+piecenum = '2'
+----setting frame, step and instance----
+frame = '1-20'
+step = '0'
+instance = '0'

source/vtkwriter_12142020.py

@@ -0,0 +1,166 @@
+INPPATH = 'C:/abq2020/test1204MC.inp'
+ODBPATH = 'C:/abq2020/test1204MC.odb'
+OUTPATH = 'C:/abq2020/test1204MC.vtk'
+
+from odbAccess import openOdb
+from abaqus import *
+from abaqusConstants import *
+
+class DataSet:
+    def __init__(self):
+        self.nodes = {}
+        self.elements = {}
+        self.stress = {}
+        self.spstress = {}
+        self.displacement = {}
+
+    def insertNode(self, nid, x, y, z):
+        self.nodes[int(nid)] = [float(x), float(y), float(z)]
+
+    def insertHexa(self, eid, n0, n1, n2, n3, n4, n5, n6, n7):
+        self.elements[int(eid)] = [int(n0), int(n1), int(n2), int(n3), int(n4), int(n5), int(n6), int(n7)]
+        
+    def insertStress(self, eid, s11, s22, s33, s12, s13, s23):
+        self.stress[int(eid)] = [float(s11), float(s22), float(s33), float(s12), float(s13), float(s23)]
+              
+    def insertSpstress(self, eid, sp1, sp2, sp3):
+        self.spstress[int(eid)] = [float(sp1), float(sp2), float(sp3)]
+
+    def insertDisplacement(self, nid, ux, uy, uz):
+        self.displacement[int(nid)] = [float(ux),  float(uy),  float(uz)]
+
+    def save(self, file_name):
+        nids = sorted(self.nodes.keys())
+        eids = sorted(self.elements.keys())
+        with open(file_name, 'w') as ofs:
+            ofs.write('<?xml version="1.0"?>')
+            ofs.write('<VTKFile type="StructuredGrid" version="0.1" byte_order="LittleEndian" header_type="UInt32">')
+            ofs.write('<StructuredGrid> WholeExtent="x1 x2 x3 y1 y2 y3 z1 z2 z3">')
+            ofs.write('<Piece Extent="x1 x2 x3 y1 y2 y3 z1 z2 z3">')
+            
+### define Points element ###
+            ofs.write('<Points>')
+            ofs.write('<DataArray type="Float32" Name="Points" NumberOfComponents="3" format="ascii">')
+            for key in nids:
+                x, y, z = self.nodes[key]
+                ofs.write("%f %f %f " % (x, y, z))                        
+            ofs.write('</DataArray>')
+            ofs.write('</Points>')  
+                      
+### define PointsData element ###            
+            ofs.write('<PointData Vectors="displacement">')
+            ofs.write('<DataArray type="Float32" NumberOfComponents="3" format="ascii" Name="displacement">')
+            for key in nids:
+                ux, uy, uz = self.displacement[key]
+                ofs.write("%f %f %f " % (ux, uy, uz))
+            ofs.write('</DataArray>')
+            ofs.write('</PointData>')
+            
+### define Cells element ###             
+            ofs.write('<Cells>')
+            ofs.write('<DataArray type="Int64" Name="connectivity" format="ascii">')
+            for key in eids:
+                ns = self.elements[key]
+                nns = map(lambda x: str(x - 1), ns)
+                cons = " ".join(nns)
+                ofs.write(cons + ' ')
+            ofs.write('</DataArray>')
+            ofs.write('<DataArray type="Int64" Name="offsets" format="ascii">')
+            for key in eids:
+                ofs.write("%d " % (key * 8))
+            ofs.write('</DataArray>')
+            ofs.write('<DataArray type="UInt8" Name="types" format="ascii">')
+            for key in eids:
+                ofs.write("10 ")
+            ofs.write('</DataArray>')
+            ofs.write('</Cells>')
+            
+### define CellData element ###              
+            ofs.write('<CellData Vectors="PrincipalStress" Tensors="Stress">')                
+            ofs.write('<DataArray type="Float32" Name="PrincipalStress" NumberOfComponents="3" format="ascii">')
+            for key in eids:
+                sp1, sp2, sp3 = self.spstress[key][0:]
+                ofs.write("%f %f %f " % (sp1, sp2, sp3))      
+            ofs.write('</DataArray>')
+            ofs.write('<DataArray type="Float32" Name="Stress" NumberOfComponents="9" format="ascii">')
+            for key in eids:
+                s11, s22, s33, s12, s13, s23 = self.stress[key][0:]
+                ofs.write("%f %f %f %f %f %f %f %f %f " % (s11, s12, s13, s12, s22, s23, s13, s23, s33))
+            ofs.write('</DataArray>')
+            ofs.write('</CellData>')
+            ofs.write('</Piece>')
+            ofs.write('</StructuredGrid>')
+            ofs.write('</VTKFile>')
+            ofs.close()
+
+
+data = DataSet()
+# 0 - nothing
+# 1 - node
+# 2 - hexahedron
+state = 0
+# read input file
+with open(INPPATH, 'r') as input_file:
+    for line in input_file:
+        if line.strip() == '*Node':
+            state = 1
+        elif line.startswith('*Element, type=C3D8R'):
+            state = 2
+        else:
+            if line.startswith("*") and (not line.startswith("**")):
+                state = 0
+            elif state == 1:
+                uid, x, y, z = line.strip().split(',')
+                data.insertNode(uid, x, y, z)
+            elif state == 2:
+                el, n0, n1, n2, n3, n4, n5, n6, n7 = line.strip().split(',')
+                data.insertHexa(el, n0, n1, n2, n3, n4, n5, n6, n7)
+            else:
+                pass
+
+# Open the odb
+myodb = openOdb(ODBPATH)
+
+# Get the frame repository for the step, find number of frames (starts at frame 0)
+stepName = myodb.steps.keys()[0]
+frames = myodb.steps[stepName].frames
+numFrames = len(frames)
+print("num Frame %d" % (numFrames))
+
+# Isolate the instance, get the number of nodes and elements
+instanceName = myodb.rootAssembly.instances.keys()[0]
+myInstance = myodb.rootAssembly.instances[instanceName]
+numNodes = len(myInstance.nodes)
+numElements = len(myInstance.elements)
+print("num Element %d, num Node %d" % (numElements, numNodes))
+
+S=myodb.steps[stepName].frames[-1].fieldOutputs['S'].getSubset(position=INTEGRATION_POINT)
+SPS=myodb.steps[stepName].frames[-1].fieldOutputs['S'].getSubset(position=INTEGRATION_POINT)
+DISP=myodb.steps[stepName].frames[-1].fieldOutputs['U']
+print("S: %d DISP: %d" % (len(S.values), len(DISP.values)))
+
+for sts in S.values:
+    eid = sts.elementLabel
+    sxx = sts.data[0]
+    syy = sts.data[1]
+    szz = sts.data[2]
+    sxy = sts.data[3]
+    sxz = sts.data[4]
+    syz = sts.data[5]
+    data.insertStress(eid, float(sxx), float(syy), float(szz), float(sxy), float(sxz), float(syz))
+
+for spsts in SPS.values:
+    eid = spsts.elementLabel  
+    smin = spsts.maxPrincipal
+    smid = spsts.midPrincipal
+    smax = spsts.minPrincipal
+    data.insertSpstress(eid, float(smin), float(smid), float(smax))
+    
+for disp in DISP.values:
+    nid = disp.nodeLabel
+    ux = disp.data[0]
+    uy = disp.data[1]
+    uz = disp.data[2]
+    data.insertDisplacement(nid, ux, uy, uz)
+
+data.save(OUTPATH)

thermal/dflux.for

@@ -0,0 +1,27 @@
+        SUBROUTINE DFLUX(FLUX,SOL,KSTEP,KINC,TIME,NOEL,NPT,COORDS,JLTYP, TEMP, PRESS, SNAME)
+        INCLUDE 'ABA_PARAM.INC'
+
+        DIMENSION COORDS(3),FLUX(2),TIME(2)
+        CHARACTER*80 SNAME
+        real U,AI,v,yita,R0,qm
+
+        U=16.5
+        AI=60.
+        v=250./60.
+        yita=0.75
+        R0=5.
+
+        qm=yita*U*AI*1000./3.14/R0/R0
+        dx=v*TIME(1)
+        dy=0.
+
+        rr= (COORDS(1)-dx)**2 + (COORDS(2)-dy)**2
+        FLUX(1)=3.*qm*exp(-3.*rr/R0/R0)
+
+        ! write(*,*) time,dx
+        ! write(*,*) COORDS,rr
+        ! write(*,*) FLUX
+        ! write(*,*)
+
+        RETURN
+        END

thermal/optx_odb2vtk.py

@@ -0,0 +1,145 @@
+from odbAccess import openOdb
+from abaqus import *
+from abaqusConstants import *
+from utils_odb2vtk import Physical_Quantity
+
+
+class CAE_Model(Physical_Quantity, object):
+    def __init__(self):
+        super(CAE_Model, self).__init__()
+        # print('self.var_type', self.var_type)
+        self.nodes = {}
+        self.elements = {}
+        self.physical_quantity = {}
+
+    def read_inp(self, INPPATH):
+        # 0 - nothing
+        # 1 - node
+        # 2 - hexahedron
+        state = 0
+        # read input file
+        with open(INPPATH, 'r') as input_file:
+            for line in input_file:
+                if line.strip() == '*Node':
+                    state = 1
+                elif line.startswith('*Element, type=C3D8R'):
+                    state = 2
+                else:
+                    if line.startswith("*") and (not line.startswith("**")):
+                        state = 0
+                    elif state == 1:
+                        nid, x, y, z = line.strip().split(',')
+                        # insertNode(uid, x, y, z)
+                        self.nodes[int(nid)] = [float(x), float(y), float(z)]
+                    elif state == 2:
+                        eid, n0, n1, n2, n3, n4, n5, n6, n7 = line.strip().split(',')
+                        # insertHexa(el, n0, n1, n2, n3, n4, n5, n6, n7)
+                        self.elements[int(eid)] = [int(n0), int(n1), int(
+                            n2), int(n3), int(n4), int(n5), int(n6), int(n7)]
+                    else:
+                        pass
+
+    def read_odb(self, ODBPATH):
+        # Open the odb
+        myodb = openOdb(ODBPATH)
+
+        # Get the frame repository for the step, find number of frames (starts at frame 0)
+        stepName = myodb.steps.keys()[0]
+        frames = myodb.steps[stepName].frames
+        numFrames = len(frames)
+        print("num Frame %d" % (numFrames))
+
+        # Isolate the instance, get the number of nodes and elements
+        instanceName = myodb.rootAssembly.instances.keys()[0]
+        myInstance = myodb.rootAssembly.instances[instanceName]
+        numNodes = len(myInstance.nodes)
+        numElements = len(myInstance.elements)
+        print("num Element %d, num Node %d" % (numElements, numNodes))
+
+        for var_id in ['U', 'A', 'V', 'RF', 'S', 'LE', 'PE', 'PEEQ', 'NT11', 'HFL', 'RFL11']:
+            # try:
+            if var_id == 'NT11':
+                # S=myodb.steps[stepName].frames[-1].fieldOutputs['S'].getSubset(position=INTEGRATION_POINT)
+                var_data = myodb.steps[stepName].frames[-1].fieldOutputs[var_id].values
+                print(var_id+": %d" % (len(var_data)))
+
+                self.pq_update(var_id)
+                self.physical_quantity[var_id] = self.insert(var_data)
+            # except:
+            #     print('jump ', var_id)
+
+    def write_vtk(self, OUTPATH, var_ids):
+        self.nids = sorted(self.nodes.keys())
+        self.eids = sorted(self.elements.keys())
+        with open(OUTPATH, 'w') as ofs:
+            ofs.write('<?xml version="1.0"?>')
+            ofs.write(
+                '<VTKFile type="StructuredGrid" version="0.1" byte_order="LittleEndian" header_type="UInt32">')
+            ofs.write(
+                '<StructuredGrid> WholeExtent="x1 x2 x3 y1 y2 y3 z1 z2 z3">')
+            ofs.write('<Piece Extent="x1 x2 x3 y1 y2 y3 z1 z2 z3">')
+
+            ### define Points element ###
+            ofs.write('<Points>')
+            ofs.write(
+                '<DataArray type="Float32" Name="Points" NumberOfComponents="3" format="ascii">')
+            for key in self.nids:
+                x, y, z = self.nodes[key]
+                ofs.write("%f %f %f " % (x, y, z))
+            ofs.write('</DataArray>')
+            ofs.write('</Points>')
+
+            ### define Cells element ###
+            ofs.write('<Cells>')
+            ofs.write(
+                '<DataArray type="Int64" Name="connectivity" format="ascii">')
+            for key in self.eids:
+                ns = self.elements[key]
+                nns = map(lambda x: str(x - 1), ns)
+                cons = " ".join(nns)
+                ofs.write(cons + ' ')
+            ofs.write('</DataArray>')
+            ofs.write('<DataArray type="Int64" Name="offsets" format="ascii">')
+            for key in self.eids:
+                ofs.write("%d " % (key * 8))
+            ofs.write('</DataArray>')
+            ofs.write('<DataArray type="UInt8" Name="types" format="ascii">')
+            for key in self.eids:
+                ofs.write("10 ")
+            ofs.write('</DataArray>')
+            ofs.write('</Cells>')
+
+            # node data/element data
+            for var_id in var_ids:
+                self.pq_update(var_id)
+
+                # var_name = var_stress.__name__
+                var_name = var_id
+                var_type = self.var_type
+                NumberOfComponents = str(self.NumberOfComponents)
+
+                ofs.write('<'+var_type+'Data '+var_type+'="' + var_name + '">')
+                ofs.write('<DataArray type="Float32" format="ascii" Name="' +
+                          var_name+'" NumberOfComponents="'+NumberOfComponents+'" >')
+                self.write(self.physical_quantity[var_id], ofs)
+                ofs.write('</DataArray>')
+                ofs.write('</'+var_type+'Data>')
+
+            # finshed
+            ofs.write('</Piece>')
+            ofs.write('</StructuredGrid>')
+            ofs.write('</VTKFile>')
+            ofs.close()
+
+
+INPPATH = 'C:/Users/Harvo/Downloads/Digital-Twin/thermal/Thermal.inp'
+ODBPATH = 'C:/Users/Harvo/Downloads/Digital-Twin/thermal/Thermal.odb'
+OUTPATH = 'C:/Users/Harvo/Downloads/Digital-Twin/thermal/Thermal.vtk'
+
+my_model = CAE_Model()
+my_model.read_inp(INPPATH)
+print('finish read_inp')
+my_model.read_odb(ODBPATH)
+print('finish read_odb')
+my_model.write_vtk(OUTPATH, var_ids=['NT11'])
+print('finish write_vtk')

thermal/utils_odb2vtk.py

@@ -0,0 +1,144 @@
+class Physical_Quantity(object):
+    def __init__(self):
+        pass
+        # Scalars NT11,RFL11
+        # Vectors displacement,HFL
+        # Tensors Stress,spstress
+
+    def pq_update(self, var_id):
+        if var_id == 'S':
+            self.pq_class = Stress()
+        if var_id == 'SPS':
+            self.pq_class = Spstress()
+        if var_id == 'DISP':
+            self.pq_class = Displacement()
+        if var_id == 'NT11':
+            self.pq_class = Temperature()
+        self.var_type = self.pq_class.var_type
+        self.NumberOfComponents = self.pq_class.NumberOfComponents
+
+    def insert(self, var_data):
+        return self.pq_class.insertData(var_data)
+
+    def write(self, *args):
+        # print(args)
+        self.pq_class.writeData(self.nids, self.eids, args)
+
+
+class Stress:
+    def __init__(self):
+        self.var_type = 'Tensors'
+        self.NumberOfComponents = 9
+
+    def insertData(self, var_data):
+        this_pq = {}
+        for var_line in var_data:
+            eid = var_line.elementLabel
+            sxx = var_line.data[0]
+            syy = var_line.data[1]
+            szz = var_line.data[2]
+            sxy = var_line.data[3]
+            sxz = var_line.data[4]
+            syz = var_line.data[5]
+            this_pq[int(eid)] = [float(sxx), float(syy), float(
+                szz), float(sxy), float(sxz), float(syz)]
+        return this_pq
+
+    def writeData(self, stress, ofs):
+        for key in self.eids:
+            s11, s22, s33, s12, s13, s23 = stress[key][0:]
+            ofs.write("%f %f %f %f %f %f %f %f %f " %
+                      (s11, s12, s13, s12, s22, s23, s13, s23, s33))
+
+
+class Spstress:
+    def __init__(self):
+        self.var_type = 'Vectors'
+        self.NumberOfComponents = 3
+
+    def insertData(self, var_data):
+        this_pq = {}
+        for var_line in var_data:
+            eid = var_line.elementLabel
+            smin = var_line.maxPrincipal
+            smin = var_line.midPrincipal
+            smin = var_line.minPrincipal
+            this_pq[int(eid)] = [float(smin), float(smin), float(smin)]
+        return this_pq
+    
+    def writeData(self, nids, eids, args):
+        spstress, ofs = args
+        for key in eids:
+            sp1, sp2, sp3 = spstress[key][0:]
+            ofs.write("%f %f %f " % (sp1, sp2, sp3))
+
+
+class Displacement:
+    def __init__(self):
+        self.var_type = 'Vectors'
+        self.NumberOfComponents = 3
+
+    def insertData(self, var_data):
+        this_pq = {}
+        for var_line in var_data:
+            nid = var_line.nodeLabel
+            ux = var_line.data[0]
+            uy = var_line.data[1]
+            uz = var_line.data[2]
+            this_pq[int(nid)] = [float(ux),  float(uy),  float(uz)]
+        return this_pq
+
+    def writeData(self, nids, eids, args):
+        displacement, ofs = args
+        for key in nids:
+            ux, uy, uz = displacement[key]
+            ofs.write("%f %f %f " % (ux, uy, uz))
+
+
+class Temperature(object):
+    def __init__(self):
+        self.var_type = 'Scalars'
+        self.NumberOfComponents = 1
+
+    def insertData(self, var_data):
+        this_pq = {}
+        for var_line in var_data:
+            nid = var_line.nodeLabel
+            temp = var_line.data
+            this_pq[int(nid)] = [float(temp)]
+        return this_pq
+
+    def writeData(self, nids, eids, args):
+        temperature, ofs = args
+        for key in nids:
+            temp = temperature[key][0]
+            # print(temp)
+            ofs.write("%f " % (temp,))
+
+
+# ### define PointsData element ###
+#             ofs.write('<PointData Vectors="displacement">')
+#             ofs.write('<DataArray type="Float32" NumberOfComponents="3" format="ascii" Name="displacement">')
+#             for key in self.nids:
+#                 ux, uy, uz = self.displacement[key]
+#                 ofs.write("%f %f %f " % (ux, uy, uz))
+#             ofs.write('</DataArray>')
+#             ofs.write('</PointData>')
+
+# ### define CellData element ###
+#             ofs.write('<CellData Vectors="PrincipalStress" Tensors="Stress">')
+#             ofs.write('<DataArray type="Float32" Name="PrincipalStress" NumberOfComponents="3" format="ascii">')
+#             for key in self.eids:
+#                 sp1, sp2, sp3 = self.spstress[key][0:]
+#                 ofs.write("%f %f %f " % (sp1, sp2, sp3))
+#             ofs.write('</DataArray>')
+#             ofs.write('<DataArray type="Float32" Name="Stress" NumberOfComponents="9" format="ascii">')
+#             for key in self.eids:
+#                 s11, s22, s33, s12, s13, s23 = self.stress[key][0:]
+#                 ofs.write("%f %f %f %f %f %f %f %f %f " % (s11, s12, s13, s12, s22, s23, s13, s23, s33))
+#             ofs.write('</DataArray>')
+#             ofs.write('</CellData>')
+#             ofs.write('</Piece>')
+#             ofs.write('</StructuredGrid>')
+#             ofs.write('</VTKFile>')
+#             ofs.close()