import sys import Vec3D FP_HAVE_NPOS = 0x00000001 FP_HAVE_KS = 0x00000002 FP_HAVE_JS = 0x00000004 FP_HAVE_TS = 0x00000008 FP_HAVE_TEMP = 0x00000010 FP_HAVE_COORD = 0x00000020 FP_HAVE_AXES = 0x00000040 FP_HAVE_KINK_ANGLE = 0x00000080 FP_HAVE_EXT = 0x00000100 FP_HAVE_NODE_ID = 0x00000200 FP_HAVE_SECT_ID = 0x00000400 FP_HAVE_GS = 0x00000800 FP_HAVE_CODS = 0x00001000 FP_HAVE_DC_POINT = 0x00002000 FP_HAVE_MODULI = 0x00004000 FP_HAVE_YS = 0x00008000 FP_HAVE_DEPTH = 0x00010000 FP_HAVE_EPJS = 0x00020000 FP_HAVE_FRAC_EXT = 0x00040000 FP_HAVE_FRAC_ANG = 0x00080000 FP_HAVE_CGM_CODE = 0x00100000 FP_HAVE_TGM_CODE = 0x00200000 CF_HAVE_ID = 0x00000001 CF_HAVE_START = 0x00000002 CF_HAVE_STOP = 0x00000004 CF_HAVE_FITFRT = 0x00000008 CF_HAVE_OLDFIT = 0x00000010 CF_HAVE_FITEXT = 0x00000020 CF_HAVE_FITPARAMS = 0x00000040 CF_HAVE_SRC_FILE = 0x00000080 CF_HAVE_EXT_MULT = 0x00000100 CS_HAVE_NAME = 0x00000001 CS_HAVE_USER_INDX = 0x00000002 CS_HAVE_LC_INFO = 0x00000004 CS_HAVE_EXT_TYPE = 0x00000008 CS_HAVE_EXT_CYCLES = 0x00000010 CS_HAVE_EXT_TIME = 0x00000020 CS_HAVE_EXT_FRAC = 0x00000040 CS_HAVE_EXT_STOP = 0x00000080 CS_HAVE_EXT_LEN = 0x00000100 CS_HAVE_TOT_PASS = 0x00000200 CS_HAVE_SIF_FLAGS = 0x00000400 CG_HAVE_START_POINT = 0x00000001 CG_HAVE_START_LEN = 0x00000002 CS_MEDIAN_EXT = 0 CS_CYCLES_EXT = 1 CS_TIME_EXT = 2 CS_USER_EXT = 3 CS_QUASI_STATIC = 4 CS_DIST_EXT = 5 STEP_LENGTH = 0 FULL_SCHED = 1 MAX_CYCLES = 2 MAX_TIME = 3 K_CRIT = 4 THRESHOLD = 5 MAX_DEPTH = 6 HCF_THRESH = 7 MAX_STEPS = 8 SIF_M_INTEGRAL = 0x00000000 SIF_DISP_CORR = 0x00000001 SIF_VCCT = 0x00000002 SIF_TYPE_MASK = 0x0000000F SIF_THERMAL_TERMS = 0x00000010 SIF_CF_TRACTIONS = 0x00000020 SIF_CF_CONTACT = 0x00000040 SIF_LARGE_ROT = 0x00000080 SIF_ELAST_PLAST = 0x00000100 SIF_INIT_STRESS = 0x00000200 SIF_INIT_STRAIN = 0x00000400 KINK_EXTEN_POLY = 0 FIXED_ORDER_POLY = 1 TRIAL_ORDER_POLY = 2 CUBIC_SPLINE = 3 MOVING_POLY = 4 NOFIT_EXTRAP = 5 NOFIT_NOEXTRAP = 6 HERMITIAN = 7 SIF_M_INTEGRAL = 0 SIF_DISP_CORR = 1 SIF_VCCT = 2 class FileIterator(object): def __init__(self,fd): self.fd = fd self.Last = None self.Next() def Next(self): if self.Last != None: self.buff = self.Last else: self.buff = self.fd.readline() if len(self.buff) > 0: self.split = self.buff.split() else: self.split = [] #print self.buff def Push(self): self.Last = self.buff def More(self): return len(self.split) > 0 def Len(self): return len(self.split) def Buff(self): return self.buff def __getitem__(self,indx): return self.split[indx] def __len__(self): return len(self.split) # ----------------------------------------------------------------------- class PointData(object): def __init__(self): self.Flags = 0 self.Ks = [] self.Gs = [] self.Js = [] self.Ts = [] self.Temps = [] self.Coord = None self.Axes = None self.KinkAngle = None self.Extension = None self.NodeId = None self.SectionId = [None,None] self.CODs = [] self.DcPoint = None self.Moduli = None self.YieldStress = None self.Depth = None self.EPJs = [] self.FracExt = [] self.FracAng = [] self.CyclesGrowthModelCode = None self.TimeGrowthModelCode = None def Read(self,fr): fr.Next() while fr.More(): if fr[0] == "NPos:": self.NPos = float(fr[1]) self.Flags |= FP_HAVE_NPOS elif fr[0] == "Ks:": num = int((fr.Len()-1)/3) for i in range(num): self.Ks.append(Vec3D.Vec3D(float(fr[i*3+1]), float(fr[i*3+2]), float(fr[i*3+3]))) self.Flags |= FP_HAVE_KS elif fr[0] == "Gs:": num = int((fr.Len()-1)/3) for i in range(num): self.Gs.append(Vec3D.Vec3D(float(fr[i*3+1]), float(fr[i*3+2]), float(fr[i*3+3]))) self.Flags |= FP_HAVE_GS elif fr[0] == "CODs:": num = int((fr.Len()-1)/3) for i in range(num): self.CODs.append(Vec3D.Vec3D(float(fr[i*3+1]), float(fr[i*3+2]), float(fr[i*3+3]))) self.Flags |= FP_HAVE_CODS elif fr[0] == "Js:": for i in range(1,fr.Len()): self.Js.append(float(fr[i])) self.Flags |= FP_HAVE_JS elif fr[0] == "Ts:": for i in range(1,fr.Len()): self.Ts.append(float(fr[i])) self.Flags |= FP_HAVE_TS elif fr[0] == "Temp:": for i in range(1,fr.Len()): self.Temps.append(float(fr[i])) self.Flags |= FP_HAVE_TEMP elif fr[0] == "Coord:": self.Coord = Vec3D.Vec3D(float(fr[1]),float(fr[2]),float(fr[3])) self.Flags |= FP_HAVE_COORD elif fr[0] == "DcPoint:": self.DcPoint = Vec3D.Vec3D(float(fr[1]),float(fr[2]),float(fr[3])) self.Flags |= FP_HAVE_DC_POINT elif fr[0] == "Axes:": self.Axes = ((Vec3D.Vec3D(float(fr[1]),float(fr[2]),float(fr[3])), Vec3D.Vec3D(float(fr[4]),float(fr[5]),float(fr[6])), Vec3D.Vec3D(float(fr[7]),float(fr[8]),float(fr[9])))) self.Flags |= FP_HAVE_AXES elif fr[0] == "KAng:": self.KinkAngle = float(fr[1]) self.Flags |= FP_HAVE_KINK_ANGLE elif fr[0] == "Ext:": self.Extension = float(fr[1]) self.Flags |= FP_HAVE_EXT elif fr[0] == "NodeId:": self.NodeId = int(fr[1]) self.Flags |= FP_HAVE_NODE_ID elif fr[0] == "SectId:": self.SectionId[0] = int(fr[1]) if fr.Len() > 2: self.SectionId[1] = int(fr[2]) self.Flags |= FP_HAVE_SECT_ID elif fr[0] == "Moduli:": self.Moduli = [float(fr[1]),float(fr[2])] self.Flags |= FP_HAVE_MODULI elif fr[0] == "YS:": self.YieldStress = float(fr[1]) self.Flags |= FP_HAVE_YS elif fr[0] == "DEPTH:": self.Depth = float(fr[1]) self.Flags |= FP_HAVE_YS elif fr[0] == "EPJs:": for i in range(1,fr.Len()): self.EPJs.append(float(fr[i])) self.Flags |= FP_HAVE_EPJS elif fr[0] == "FExt:": for i in range(1,fr.Len()): self.FracExt.append(fr[i]) self.Flags |= FP_HAVE_FRAC_EXT elif fr[0] == "FAng:": for i in range(1,fr.Len()): self.FracAng.append(fr[i]) self.Flags |= FP_HAVE_FRAC_ANG elif fr[0] == "CGM:": self.CyclesGrowthModelCode = fr[1] self.Flags |= FP_HAVE_CGM_CODE elif fr[0] == "TGM:": self.TimeGrowthModelCode = fr[1] self.Flags |= FP_HAVE_TGM_CODE elif fr[0] == "ENDFP": break fr.Next() def Save(self,out): if (self.Flags & FP_HAVE_NPOS) != 0: print(" NPos:",self.NPos, file=out) if (self.Flags & FP_HAVE_KS) != 0: print(" Ks:", end=' ', file=out) for k in self.Ks: print(k[0],k[1],k[2], end=' ', file=out) print(file=out) if (self.Flags & FP_HAVE_GS) != 0: print(" Gs:", end=' ', file=out) for g in self.Gs: print(g[0],g[1],g[2], end=' ', file=out) print(file=out) if (self.Flags & FP_HAVE_CODS) != 0: print(" CODs:", end=' ', file=out) for c in self.CODs: print(c[0],c[1],c[2], end=' ', file=out) print(file=out) if (self.Flags & FP_HAVE_JS) != 0: print(" Js:", end=' ', file=out) for j in self.Js: print(j, end=' ', file=out) print(file=out) if (self.Flags & FP_HAVE_TS) != 0: print(" Ts:", end=' ', file=out) for t in self.Ts: print(t, end=' ', file=out) print(file=out) if (self.Flags & FP_HAVE_TEMP) != 0: print(" Temp:", end=' ', file=out) for t in self.Temps: print(t, end=' ', file=out) print(file=out) if (self.Flags & FP_HAVE_COORD) != 0: print(" Coord:",self.Coord[0],self.Coord[1],self.Coord[2], file=out) if (self.Flags & FP_HAVE_DC_POINT) != 0: print(" DcPoint:",self.DcPoint[0],self.DcPoint[1],self.DcPoint[2], file=out) if (self.Flags & FP_HAVE_AXES) != 0: print(" Axes:", \ self.Axes[0][0],self.Axes[0][1],self.Axes[0][2], \ self.Axes[1][0],self.Axes[1][1],self.Axes[1][2], \ self.Axes[2][0],self.Axes[2][1],self.Axes[2][2], file=out) if (self.Flags & FP_HAVE_KINK_ANGLE) != 0: print(" KAng:",self.KinkAngle, file=out) if (self.Flags & FP_HAVE_EXT) != 0: print(" Ext:",self.Extension, file=out) if (self.Flags & FP_HAVE_NODE_ID) != 0: print(" NodeId:",self.NodeId, file=out) if (self.Flags & FP_HAVE_SECT_ID) != 0: if self.SectionId[1]: print(" SectId:",self.SectionId[0],self.SectionId[1], file=out) else: print(" SectId:",self.SectionId[0], file=out) if (self.Flags & FP_HAVE_MODULI) != 0: print(" Moduli:",self.Moduli[0],self.Moduli[1], file=out) if (self.Flags & FP_HAVE_YS) != 0: print(" YS:",self.YieldStress, file=out) if (self.Flags & FP_HAVE_DEPTH) != 0: print(" DEPTH:",self.Depth, file=out) if (self.Flags & FP_HAVE_EPJS) != 0: print(" EPJs:", end=' ', file=out) for j in self.EPJs: print(j, end=' ', file=out) print(file=out) if (self.Flags & FP_HAVE_FRAC_EXT) != 0: print(" FExt:", end=' ', file=out) for j in self.FracExt: print(j, end=' ', file=out) print(file=out) if (self.Flags & FP_HAVE_FRAC_ANG) != 0: print(" FAng:", end=' ', file=out) for j in self.FracAng: print(j, end=' ', file=out) print(file=out) if (self.Flags & FP_HAVE_CGM_CODE) != 0: print(" CGM:",self.CyclesGrowthModelCode, file=out) if (self.Flags & FP_HAVE_TGM_CODE) != 0: print(" CGM:",self.TimeGrowthModelCode, file=out) print("ENDFP", file=out) # ----------------------------------------------------------------------- class FitParams(object): def __init__(self): self.SmoothType = None self.PolyOrder = None self.Discard = [None,None] self.Extrapolate = [None,None] class FrontData(object): def __init__(self): self.Flags = 0 self.Id = None self.StartPoint = None self.StopPoint = None self.FitFront = [] self.FitOldFront = [] self.FitExtens = [] self.FitParams = FitParams() self.SrcFile = None self.ExtensionMult = None self.PointData = [] def Read(self,fr): fr.Next() while fr.More(): if fr[0] == "Id:": self.Id = int(fr[1]) self.Flags |= CF_HAVE_ID elif fr[0] == "Start:": self.StartPoint = Vec3D.Vec3D(float(fr[1]),float(fr[2]),float(fr[3])) self.Flags |= CF_HAVE_START elif fr[0] == "Stop:": self.StopPoint = Vec3D.Vec3D(float(fr[1]),float(fr[2]),float(fr[3])) self.Flags |= CF_HAVE_STOP elif fr[0] == "FitFrt:": num = int((fr.Len()-1)/3) for i in range(num): self.FitFront.append(Vec3D.Vec3D(float(fr[i*3+1]), float(fr[i*3+2]), float(fr[i*3+3]))) self.Flags |= CF_HAVE_FITFRT elif fr[0] == "FitOld:": num = int((fr.Len()-1)/3) for i in range(num): self.FitOldFront.append(Vec3D.Vec3D(float(fr[i*3+1]), float(fr[i*3+2]), float(fr[i*3+3]))) self.Flags |= CF_HAVE_OLDFIT elif fr[0] == "FitExt:": for i in range(1,fr.Len()): self.FitExtens.append(float(fr[i])) self.Flags |= CF_HAVE_FITEXT elif fr[0] == "FitParams:": if fr[1] == "KINK_EXTEN_POLY": self.FitParams.SmoothType = KINK_EXTEN_POLY elif fr[1] == "FIXED_ORDER_POLY": self.FitParams.SmoothType = FIXED_ORDER_POLY elif fr[1] == "MULTIPLE_POLY": self.FitParams.SmoothType = MULTIPLE_POLY elif fr[1] == "CUBIC_SPLINE": self.FitParams.SmoothType = CUBIC_SPLINE elif fr[1] == "MOVING_POLY": self.FitParams.SmoothType = MOVING_POLY elif fr[1] == "NOFIT_EXTRAP": self.FitParams.SmoothType = NOFIT_EXTRAP elif fr[1] == "NOFIT_NOEXTRAP": self.FitParams.SmoothType = NOFIT_NOEXTRAP elif fr[1] == "HERMITIAN": self.FitParams.SmoothType = HERMITIAN elif fr[1] == "PARTIAL": self.FitParams.SmoothType = PARTIAL self.FitParams.PolyOrder = int(fr[2]) self.FitParams.Discard[0] = int(fr[3]) self.FitParams.Discard[1] = int(fr[4]) self.FitParams.Extrapolate[0] = float(fr[5]) self.FitParams.Extrapolate[1] = float(fr[6]) self.Flags |= CF_HAVE_FITPARAMS elif fr[0] == "SrcFile:": self.SrcFile = fr[1] for i in range(2,fr.Len()): self.SrcFile += ' ' self.SrcFile += fr[i] self.Flags |= CF_HAVE_SRC_FILE elif fr[0] == "FrtPts:": num = int(fr[1]) for i in range(num): self.PointData.append(PointData()) self.PointData[-1].Read(fr) elif fr[0] == "ExtMult:": self.ExtensionMult = float(fr[1]) self.Flags |= CF_HAVE_EXT_MULT elif fr[0] == "ENDCF": break fr.Next() def Save(self,out): if (self.Flags & CF_HAVE_ID) != 0: print(" ID:",self.Id, file=out) if (self.Flags & CF_HAVE_START) != 0: print(" Start:",self. \ StartPoint[0],self.StartPoint[1],self.StartPoint[2], file=out) if (self.Flags & CF_HAVE_STOP) != 0: print(" Stop:", \ self.StopPoint[0],self.StopPoint[1],self.StopPoint[2], file=out) if (self.Flags & CF_HAVE_FITFRT) != 0: print(" FitFrt: ", end=' ', file=out) for i in range(len(self.FitFront)): print(self.FitFront[i][0], \ self.FitFront[i][1], \ self.FitFront[i][2], end=' ', file=out) print(file=out) if (self.Flags & CF_HAVE_OLDFIT) != 0: print(" FitOld:", end=' ', file=out) for i in range(len(self.FitOldFront)): print(self.FitOldFront[i][0], \ self.FitOldFront[i][1], \ self.FitOldFront[i][2], end=' ', file=out) print(file=out) if (self.Flags & CF_HAVE_FITEXT) != 0: print(" FitExt:", end=' ', file=out) for i in range(len(self.FitOldFront)): print(self.FitExtens[i], end=' ', file=out) print(file=out) if (self.Flags & CF_HAVE_FITPARAMS) != 0: print(" FitParams:", end=' ', file=out) if self.FitParams.SmoothType == KINK_EXTEN_POLY: print("KINK_EXTEN_POLY", end=' ', file=out) elif self.FitParams.SmoothType == FIXED_ORDER_POLY: print("FIXED_ORDER_POLY", end=' ', file=out) elif self.FitParams.SmoothType == MULTIPLE_POLY: print("MULTIPLE_POLY", end=' ', file=out) elif self.FitParams.SmoothType == CUBIC_SPLINE: print("CUBIC_SPLINE", end=' ', file=out) elif self.FitParams.SmoothType == MOVING_POLY: print("MOVING_POLY", end=' ', file=out) elif self.FitParams.SmoothType == NOFIT_EXTRAP: print("NOFIT_EXTRAP", end=' ', file=out) elif self.FitParams.SmoothType == NOFIT_NOEXTRAP: print("NOFIT_NOEXTRAP", end=' ', file=out) elif self.FitParams.SmoothType == HERMITIAN: print("HERMITIAN", end=' ', file=out) elif self.FitParams.SmoothType == PARTIAL: print("PARTIAL", end=' ', file=out) print(self.FitParams.PolyOrder,self.FitParams.Discard[0], \ self.FitParams.Discard[1],self.FitParams.Extrapolate[0], \ self.FitParams.Extrapolate[1], file=out) if (self.Flags & CF_HAVE_SRC_FILE) != 0: print(" SrcFile:",self.SrcFile, file=out) if (self.Flags & CF_HAVE_EXT_MULT) != 0: print(" ExtMult: ",self.ExtensionMult, file=out) print(" FrtPts:",len(self.PointData), file=out) for point in self.PointData: point.Save(out) print("ENDCF", file=out) # ----------------------------------------------------------------------- class Step(object): def __init__(self): self.Flags = 0 self.Name = None self.UserIndx = None self.ExtensionType = None self.ExtensionCycles = None self.ExtensionTime = None self.ExtensionFraction = None self.ExtensionLength = None self.TotalPasses = None self.StopReason = None self.SifComputationFlags = None self.FrontData = [] def Read(self,fr): fr.Next() while fr.More(): if fr[0] == "Name:": self.Name = fr[1] for i in range(2,fr.Len()): self.Name += ' ' self.Name += fr[i] self.Flags |= CS_HAVE_NAME elif fr[0] == "Extension:": if fr[1] == "MEDIAN": self.ExtensionType = CS_MEDIAN_EXT elif fr[1] == "DIST": self.ExtensionType = CS_DIST_EXT elif fr[1] == "CYCLES": self.ExtensionType = CS_CYCLES_EXT elif fr[1] == "TIME": self.ExtensionType = CS_TIME_EXT self.Flags |= CS_HAVE_EXT_TYPE elif fr[0] == "ExType:": if fr[1] == "MEDIAN": self.ExtensionType = CS_MEDIAN_EXT elif fr[1] == "DIST": self.ExtensionType = CS_DIST_EXT elif fr[1] == "CYCLES": self.ExtensionType = CS_CYCLES_EXT elif fr[1] == "TIME": self.ExtensionType = CS_TIME_EXT elif fr[1] == "USER": self.ExtensionType = CS_USER_EXT elif fr[1] == "QUASI_STATIC": self.ExtensionType = CS_QUASI_STATIC self.Flags |= CS_HAVE_EXT_TYPE elif fr[0] == "Cycles:": self.ExtensionCycles = float(fr[1]) self.Flags |= CS_HAVE_EXT_CYCLES elif fr[0] == "Time:": self.ExtensionTime = float(fr[1]) self.Flags |= CS_HAVE_EXT_TIME elif fr[0] == "Fraction:": self.ExtensionFraction = float(fr[1]) self.Flags |= CS_HAVE_EXT_FRAC elif fr[0] == "Length:": self.ExtensionLength = float(fr[1]) self.Flags |= CS_HAVE_EXT_LEN elif fr[0] == "TotalPasses:": self.TotalPasses = float(fr[1]) self.Flags |= CS_HAVE_TOT_PASS elif fr[0] == "StopReason:": if fr[1] == "STEP_LENGTH": self.StopReason = STEP_LENGTH elif fr[1] == "FULL_SCHED": self.StopReason = FULL_SCHED elif fr[1] == "MAX_CYCLES": self.StopReason = MAX_CYCLES elif fr[1] == "MAX_TIME": self.StopReason = MAX_TIME elif fr[1] == "K_CRIT": self.StopReason = K_CRIT elif fr[1] == "THRESHOLD": self.StopReason = THRESHOLD elif fr[1] == "MAX_DEPTH": self.StopReason = MAX_DEPTH elif fr[1] == "HCF_THRESH": self.StopReason = HCF_THRESH self.Flags |= CS_HAVE_EXT_STOP elif fr[0] == "SifCompType:": if fr[1] == "M_INTEGRAL": self.SifComputationType = SIF_M_INTEGRAL elif fr[1] == "DISP_CORR": self.SifComputationType = SIF_DISP_CORR elif fr[1] == "VCCT": self.SifComputationType = SIF_VCCT self.Flags |= CS_HAVE_SIF_FLAGS elif fr[0] == "SifCompFlags:": if fr[1] == "M_INTEGRAL": self.SifComputationFlags = SIF_M_INTEGRAL elif fr[1] == "DISP_CORR": self.SifComputationFlags = SIF_DISP_CORR elif fr[1] == "VCCT": self.SifComputationFlags = SIF_VCCT for i in range(2,fr.Len()): if fr[i] == "THERM": self.SifComputationFlags |= SIF_THERMAL_TERMS elif fr[i] == "CFT": self.SifComputationFlags |= SIF_CF_TRACTIONS elif fr[i] == "CFC": self.SifComputationFlags |= SIF_CF_CONTACT elif fr[i] == "LROT": self.SifComputationFlags |= SIF_LARGE_ROT elif fr[i] == "EPJ": self.SifComputationFlags |= SIF_ELAST_PLAST elif fr[i] == "SIG0": self.SifComputationFlags |= SIF_INIT_STRESS elif fr[i] == "EPS0": self.SifComputationFlags |= SIF_INIT_STRAIN self.Flags |= CS_HAVE_SIF_FLAGS elif fr[0] == "LoadFactors:": num = int(fr[1]) for i in range(num): fr.Next() elif fr[0] == "FrtData:": num = int(fr[1]) for i in range(num): self.FrontData.append(FrontData()) self.FrontData[-1].Read(fr) elif fr[0] == "ENDCS": break fr.Next() def Save(self,out): if (self.Flags & CS_HAVE_NAME) != 0: print(" Name:",Name, file=out) if (self.Flags & CS_HAVE_EXT_TYPE) != 0: print(" ExtType:", end=' ', file=out) if self.ExtensionType == CS_MEDIAN_EXT: print("MEDIAN", file=out) elif self.ExtensionType == CS_DIST_EXT: print("DIST",file=out) elif self.ExtensionType == CS_CYCLES_EXT: print("CYCLES", file=out) elif self.ExtensionType == CS_TIME_EXT: print("TIME", file=out) elif self.ExtensionType == CS_USER_EXT: print("USER", file=out) elif self.ExtensionType == CS_QUASI_STATIC: print("QUASI_STATIC", file=out) if (self.Flags & CS_HAVE_EXT_CYCLES) != 0: print(" Cycles:",self.ExtensionCycles, file=out) if (self.Flags & CS_HAVE_EXT_TIME) != 0: print(" Time:",self.ExtensionTime, file=out) if (self.Flags & CS_HAVE_EXT_FRAC) != 0: print(" Fraction:",self.ExtensionFraction, file=out) if (self.Flags & CS_HAVE_EXT_LEN) != 0: print(" Length:",self.ExtensionLength, file=out) if (self.Flags & CS_HAVE_EXT_STOP) != 0: print(" StopReason:", end=' ', file=out) if self.StopReason == STEP_LENGTH: print("STEP_LENGTH", file=out) elif self.StopReason == FULL_SCHED: print("FULL_SCHED", file=out) elif self.StopReason == MAX_CYCLES: print("MAX_CYCLES", file=out) elif self.StopReason == MAX_TIME: print("MAX_TIME", file=out) elif self.StopReason == K_CRIT: print("K_CRIT", file=out) elif self.StopReason == THRESHOLD: print("THRESHOLD", file=out) elif self.StopReason == MAX_DEPTH: print("MAX_DEPTH", file=out) elif self.StopReason == HCF_THRESH: print("HCF_THRESH", file=out) if (self.Flags & CS_HAVE_SIF_FLAGS) != 0: print(" SifCompFlags:", end=' ', file=out) comp_type = self.SifComputationFlags & SIF_TYPE_MASK if comp_type == SIF_M_INTEGRAL: print("M_INTEGRAL", end=' ', file=out) elif comp_type == SIF_DISP_CORR: print("DISP_CORR", end=' ', file=out) elif comp_type == SIF_VCCT: print("VCCT", end=' ', file=out) if (self.SifComputationFlags & SIF_THERMAL_TERMS) != 0: print("THERM", end=' ', file=out) if (self.SifComputationFlags & SIF_CF_TRACTIONS) != 0: print("CFT", end=' ', file=out) if (self.SifComputationFlags & SIF_CF_CONTACT) != 0: print("CFC", end=' ', file=out) if (self.SifComputationFlags & SIF_LARGE_ROT) != 0: print("LROT", end=' ', file=out) if (self.SifComputationFlags & SIF_ELAST_PLAST) != 0: print("EPJ", end=' ', file=out) if (self.SifComputationFlags & SIF_INIT_STRESS) != 0: print("SIG0", end=' ', file=out) if (self.SifComputationFlags & SIF_INIT_STRAIN) != 0: print("EPS0", end=' ', file=out) print(file=out) print(" FrtData: ",len(self.FrontData), file=out) for front in self.FrontData: front.Save(out) ; print("ENDCS", file=out) # ----------------------------------------------------------------------- class LoadStepMap(object): class LoadStepCFTInfo(object): def __init__(self,int_step,ext_step,sub,anal_type,mesh,stress): self.IntStep = int_step self.ExtStep = ext_step self.Substep = sub self.AnalType = anal_type self.MeshFile = mesh self.StressFile = stress def __init__(self): self.Map = [] self.RevMap = [] #self.MaxSub = [] self.SubRange = [] self.Flags = [] self.Labels = {} self.ExtData = {} def Read(self,fr): fr.Next() if fr[0] != "LOAD_STEP_MAP": self._FormatError(fr) fr.Next() if fr[0] != "(": self._FormatError(fr) fr.Next() if fr[0] != "VERSION:" or len(fr) < 2: self._FormatError(fr) if int(fr[1]) > 1: self._VersionError(fr) fr.Next() if fr[0] == "SUB_RANGE:": if len(fr) < 2: self._FormatError(fr) num = int(fr[1]) fr.Next() if len(fr) < 3*num: self._FormatError(fr) for i in range(num): self.SubRange.append((int(fr[i*3]),int(fr[i*3+1]))) self.Flags.append(int(fr[i*3+2])) elif fr[0] != "MAX_SUB:": if len(fr) < 2: self._FormatError(fr) num = int(fr[1]) fr.Next() if len(fr) < 2*num: self._FormatError(fr) for i in range(num): self.SubRange.append((0,int(fr[i*2]))) self.Flags.append(int(fr[i*2+1])) else: self._FormatError(fr) fr.Next() while fr[0] != ")": if fr[0] == "LABELS:": if len(fr) < 2: self._FormatError(fr) lnum = int(fr[1]) for i in range(lnum): fr.Next() if len(fr) < 2: self._FormatError(fr) step_id = int(fr[0]) lab = fr[1] for j in range(2,len(fr)): lab += ' ' lab += fr[j] self.Labels[step_id] = lab elif fr[0] == "EXT_CFT:": if len(fr) < 2: self._FormatError(fr) cnum = int(fr[1]) for i in range(cnum): fr.Next() if len(fr) < 6: self._FormatError(fr) step_id = int(fr[0]) self.ExtData[step_id] = \ self.LoadStepCFTInfo(-1,int(fr[4]),int(fr[5]), fr[1],fr[2],fr[3]) else: self._FormatError(fr) fr.Next() def Save(self,out,indent): print(indent,"LOAD_STEP_MAP", file=out) print(indent,"(", file=out) print(indent,"VERSION: ",1, file=out) #print(indent,"MAX_SUB: ",len(self.MaxSub), file=out) #if len(self.MaxSub) > 0: # print(indent, end=' ', file=out) # for i in range(len(self.MaxSub)): # print(self.MaxSub[i],self.Flags[i], end=' ', file=out) # print(file=out) print(indent,"SUB_RANGE:", len(self.SubRange), file=out) if len(self.SubRange) > 0: print(indent, end=' ', file=out) for i in range(len(self.SubRange)): print(self.SubRange[i][0],self.SubRange[i][1], self.Flags[i], end=' ', file=out) print(file=out) if len(self.Labels) > 0: print(indent,"LABELS:",len(self.Labels), file=out) for step_id in self.Labels: print(indent,step_id,self.Labels[step_id], file=out) if len(self.ExtData) > 0: print(indent,"EXT_CFT:",len(self.ExtData), file=out) for step_id in self.ExtData: cft_info = self.ExtData[step_id] print(indent,step_id,cft_info.AnalType, \ cft_info.MeshFile,cft_info.StressFile, file=out) print(indent,")", file=out) def _FormatError(self,fr): print("Format Error:",fr.Buff(), file=sys.stderr) raise SystemExit def _VersionError(self,fr): print("Version Error:",fr.Buff(), file=sys.stderr) raise SystemExit # ----------------------------------------------------------------------- class GrowthData(object): def __init__(self,file_name=None): self.Flags = 0 self.Steps = [] self.StartPoint = None self.StartLength = None self.LoadCaseMaps = [] if file_name: self.ReadFile(file_name) def ReadFile(self,file_name): fd = open(file_name) self.ReadFileFd(fd) def ReadFileFd(self,fd): fr = FileIterator(fd) self.ReadFileFr(fr) def ReadFileFr(self,fr): self.Flags = 0 self.Steps = [] while fr.More(): if fr[0] == "StartPoint:": self.StartPoint = Vec3D.Vec3D(float(fr[1]),float(fr[2]),float(fr[3])) self.Flags |= CG_HAVE_START_POINT elif fr[0] == "StartLen:": self.StartLength = float(fr[1]) Flags |= CG_HAVE_START_LEN elif fr[0] == "StepData:": num = int(fr[1]) for i in range(num): self.Steps.append(Step()) self.Steps[-1].Read(fr) elif fr[0] == "LoadCaseMaps:": num = int(fr[1]) for i in range(num): fr.Next() step = int(fr[1]) if len(self.LoadCaseMaps) == 0 or \ step != self.LoadCaseMaps[-1][0]: self.LoadCaseMaps.append((step,LoadStepMap())) self.LoadCaseMaps[-1][1].Read(fr) ; elif fr[0] == "ENDFCGD": break fr.Next() def WriteFile(self,out,start=0,stop=-1): if (self.Flags & CG_HAVE_START_POINT) != 0: print("StartPoint:",self.StartPoint, file=out) elif (self.Flags & CG_HAVE_START_LEN) != 0: print("StartLen:",self.StartLength, file=out) last = stop if last < 0: last = len(self.Steps)-1 num = last - start + 1 #print >> sys.stderr,"first,last,num:",start,stop,num print("StepData:",num, file=out) for i in range(start,last+1): self.Steps[i].Save(out) ; print("LoadCaseMaps:",len(self.LoadCaseMaps), file=out) for i in range(len(self.LoadCaseMaps)): print(" LoadMap:",self.LoadCaseMaps[i][0], file=out) self.LoadCaseMaps[i][1].Save(out," ") print("ENDFCGD", file=out) def NumSteps(self): return len(self.Steps) def NumFronts(self,step): return len(self.Steps[step].FrontData) def NumFrontPoints(self,step,front): if (step < len(self.Steps)) and (front < len(self.Steps[step].FrontData)): return len(self.Steps[step].FrontData[front].PointData) else: return 0 # ---- Front Data def GetFitFront(self,step,front): return self.Steps[step].FrontData[front].FitFront def GetFitOldFront(self,step,front): return self.Steps[step].FrontData[front].FitOldFront # ---- Point Data def GetFrontStart(self,step,front): return self.Steps[step].FrontData[front].StartPoint def GetFrontStop(self,step,front): return self.Steps[step].FrontData[front].StopPoint def GetFrontPointNPos(self,step,front,index): return self.Steps[step].FrontData[front].PointData[index].NPos def GetFrontPointCoord(self,step,front,index): return self.Steps[step].FrontData[front].PointData[index].Coord def GetFrontPointExtension(self,step,front,index): return self.Steps[step].FrontData[front].PointData[index].Extension def GetFrontPointKinkAngleRad(self,step,front,index): return self.Steps[step].FrontData[front].PointData[index].KinkAngle def GetFrontPointAxes(self,step,front,index): return self.Steps[step].FrontData[front].PointData[index].Axes def GetFrontPointKs(self,step,front,index): return self.Steps[step].FrontData[front].PointData[index].Ks def GetFrontPointSumKs(self,step,front,index): tmp = Vec3D.Vec3D(0,0,0) for K in self.Steps[step].FrontData[front].PointData[index].Ks: tmp += K return K def SetExtensionCycles(self,step,cycles): self.Steps[step].ExtensionCycles = cycles self.Steps[step].Flags |= CS_HAVE_EXT_CYCLES # ================================================================= def GetFile(pos): if len(sys.argv) <= pos: print("No file specified!", file=sys.stderr) raise SystemExit return GrowthData(sys.argv[pos]) def FindKAtPos(pos,point_data): for i in range(len(point_data)-1): if pos >= point_data[i].NPos and pos < point_data[i+1].NPos: par = (pos - point_data[i].NPos) / \ (point_data[i+1].NPos - point_data[i].NPos) Ks = (1.0-par)*point_data[i].Ks[0] + \ par*point_data[i].Ks[0] crd = (1.0-par)*point_data[i].Coord + \ par*point_data[i].Coord return crd,Ks # ================================================================== class SifEntry: def __init__(self): self.NormCoord = None self.CartCoord = None self.CrkFrame = None self.Ks = None self.Js = None self.Ts = None self.Temp = None self.FrontIndx = None class SifHistCalc : def __init__(self,growth_data) : self.GrowthData = growth_data def FindContainingSeg(self,par,step,front_indx): p0 = self.GrowthData.GetFrontPointNPos(step,front_indx,0) p1 = self.GrowthData.GetFrontPointNPos(step,front_indx,-1) if par < p0: return 0 if par > p1: return self.GrowthData.NumFrontPoints(step,front_indx)-2 for i in range(self.GrowthData.NumFrontPoints(step,front_indx)-1): p1 = self.GrowthData.GetFrontPointNPos(step,front_indx,i+1) if par >= p0 and par <= p1: return i p0 = p1 assert(False) return None def ClosestPath(self,ndist_start,fronts): path = [] # find the start point for the path seg = self.FindContainingSeg(ndist_start,fronts[0][0],fronts[0][1]) ; p0 = self.GrowthData.GetFrontPointNPos(fronts[0][0],fronts[0][1],seg) p1 = self.GrowthData.GetFrontPointNPos(fronts[0][0],fronts[0][1],seg+1) # linear interpolate to find starting coordinates and K's par = (ndist_start-p0) / (p1-p0) self.AppendInterpolateEntry(par,fronts[0][0],fronts[0][1], seg,seg+1,path) ; # now loop to find the closest point on each subsequent crack front ndist = ndist_start ; curr_point = path[-1].CartCoord ; for i in range(1,len(fronts)): step = fronts[i][0] findx = fronts[i][1] if self.GrowthData.NumFrontPoints(step,findx) == 0: continue ; # loop through the segments min_dist = 0 min_par = 0 min_seg = 0 for i in range(self.GrowthData.NumFrontPoints(step,findx)-1): p0 = self.GrowthData.GetFrontPointCoord(step,findx,i) p1 = self.GrowthData.GetFrontPointCoord(step,findx,i+1) delta = p1-p0 loc_par = ((curr_point-p0) * delta) / (delta * delta) close = p0 + loc_par * delta if loc_par >= 0.0 and loc_par <= 1.0: dist = (close-curr_point).Magnitude() elif loc_par < 0.0: dist = (curr_point-p0).Magnitude() loc_par = 0.0 else: dist = (curr_point-p1).Magnitude() loc_par = 1.0 if i==0 or dist < min_dist: min_dist = dist min_seg = i min_par = loc_par self.AppendInterpolateEntry(min_par,step,findx, min_seg,min_seg+1,path) ; ndist = min_dist curr_point = path[-1].CartCoord ; return path def AppendInterpolateEntry(self,par,step,front,pt0,pt1,path): path.append(self.InterpolateEntry(par,step,front,pt0,pt1)) def InterpolateEntry(self,par,step,front,pt0,pt1): entry = SifEntry() entry.NormCoord = (1-par)*self.GrowthData.GetFrontPointNPos(step,front,pt0) +\ par*self.GrowthData.GetFrontPointNPos(step,front,pt1) ; entry.CartCoord = (1-par)*self.GrowthData.GetFrontPointCoord(step,front,pt0) +\ par*self.GrowthData.GetFrontPointCoord(step,front,pt1) ; b0 = self.GrowthData.GetFrontPointAxes(step,front,pt0) b1 = self.GrowthData.GetFrontPointAxes(step,front,pt1) e0 = ((1-par)*b0[0] + par*b1[0]).Normalize() e1 = ((1-par)*b0[1] + par*b1[1]).Normalize() e2 = Vec3D.CrossProd(e0,e1) entry.CrkFrame = (e0,e1,e2) ks0 = self.GrowthData.GetFrontPointKs(step,front,pt0) ; ks1 = self.GrowthData.GetFrontPointKs(step,front,pt1) ; entry.Ks = [] for i in range(len(ks0)): k = (1-par)*ks0[i] + par*ks1[i] entry.Ks.append(k) # J0 = self.GrowthData.GetFrontPointJs(step,front,pt0) # J1 = self.GrowthData.GetFrontPointJs(step,front,pt1) # entry.Js = [] # for i in range(len(J0)): # J = (1-par)*J0[i] + par*J1[i] # entry.Js.append(J) # # T0 = self.GrowthData.GetFrontPointTs(step,front,pt0) # T1 = self.GrowthData.GetFrontPointTs(step,front,pt1) # entry.Ts = [] # for i in range(len(T0)): # T = (1-par)*T0[i] + par*T1[i] # entry.Ts.append(T) # # t0 = self.GrowthData.GetFrontPointTemp(step,front,pt0) # t1 = self.GrowthData.GetFrontPointTemp(step,front,pt1) # entry.Temp = [] # for i in range(len(t0)): # temp = (1-par)*t0[i] + par*t1[i] # entry.Temp.append(temp) entry.FrontIndx = front return entry # ----------------------------------------------------------------- if __name__ == "__main__": if len(sys.argv) < 2: print("%s -step prints a list of steps in the file" % (sys.argv[0],)) print("%s -info prints a list of the data available for each step" % (sys.argv[0],)) print("%s -extract [-start x] [-stop y] extracts the steps from x to y" % (sys.argv[0],)) print("%s -starts print the front start points" % (sys.argv[0],)) print("%s -merge ... merge two or more files" % (sys.argv[0],)) print("%s -path npos prints the K's and coordinates at the normalize position" % (sys.argv[0],)) print("%s -nearest_path npos prints a and K along a nearest path" % (sys.argv[0],)) raise SystemExit if sys.argv[1] == "-step": data = GetFile(2) num = data.NumSteps() print("Steps: 0 - %d" % (num-1,)) elif sys.argv[1] == "-info": data = GetFile(2) print("Global Data:") if data.Flags & CG_HAVE_START_POINT: print(" Start point") if data.Flags & CG_HAVE_START_LEN: print(" Start length") print() print("Per Step Data:") if data.Steps[0].Flags & CS_HAVE_NAME: print(" step name") if data.Steps[0].Flags & CS_HAVE_USER_INDX: print(" user index") if data.Steps[0].Flags & CS_HAVE_LC_INFO: print(" load case info") if data.Steps[0].Flags & CS_HAVE_EXT_TYPE: print(" extension type") if data.Steps[0].Flags & CS_HAVE_EXT_CYCLES: print(" cycles") if data.Steps[0].Flags & CS_HAVE_EXT_TIME: print(" time") if data.Steps[0].Flags & CS_HAVE_EXT_FRAC: print(" growth fraction") if data.Steps[0].Flags & CS_HAVE_EXT_STOP: print(" stop reason") if data.Steps[0].Flags & CS_HAVE_EXT_LEN: print(" growth length") if data.Steps[0].Flags & CS_HAVE_SIF_FLAGS: print(" SIF computation flags") print() print("Per Front Data:") if data.Steps[0].FrontData[0].Flags & CF_HAVE_ID: print(" front id") if data.Steps[0].FrontData[0].Flags & CF_HAVE_START: print(" surface start point") if data.Steps[0].FrontData[0].Flags & CF_HAVE_STOP: print(" surface stop point") if data.Steps[0].FrontData[0].Flags & CF_HAVE_FITFRT: print(" fit points") if data.Steps[0].FrontData[0].Flags & CF_HAVE_OLDFIT: print(" unfit points") if data.Steps[0].FrontData[0].Flags & CF_HAVE_FITEXT: print(" fit extensions") if data.Steps[0].FrontData[0].Flags & CF_HAVE_FITPARAMS: print(" fit parameters") if data.Steps[0].FrontData[0].Flags & CF_HAVE_SRC_FILE: print(" source file") if data.Steps[0].FrontData[0].Flags & CF_HAVE_EXT_MULT: print(" extension multiplier") print() print("Per Point Data:") tmp = data.Steps[0].FrontData[0].PointData[0] if tmp.Flags & FP_HAVE_NPOS: print(" normalized position") if tmp.Flags & FP_HAVE_KS: print(" SIF's") if tmp.Flags & FP_HAVE_JS: print(" J-Integral") if tmp.Flags & FP_HAVE_TS: print(" T-stress") if tmp.Flags & FP_HAVE_TEMP: print(" temperatures") if tmp.Flags & FP_HAVE_COORD: print(" coordinate") if tmp.Flags & FP_HAVE_AXES: print(" local axis") if tmp.Flags & FP_HAVE_KINK_ANGLE: print(" kink angle") if tmp.Flags & FP_HAVE_EXT: print(" extension") if tmp.Flags & FP_HAVE_NODE_ID: print(" front node id") if tmp.Flags & FP_HAVE_SECT_ID: print(" front element section id") if tmp.Flags & FP_HAVE_GS: print(" Gs") if tmp.Flags & FP_HAVE_CODS: print(" CODs") if tmp.Flags & FP_HAVE_DC_POINT: print(" DC point coord") if tmp.Flags & FP_HAVE_MODULI: print(" E and nu") if tmp.Flags & FP_HAVE_YS: print(" Yield stress") if tmp.Flags & FP_HAVE_DEPTH: print(" depth") elif sys.argv[1] == "-extract": start = -1 stop = -1 for i in range(2,len(sys.argv)): if sys.argv[i] == "-start": start = int(sys.argv[i+1]) elif sys.argv[i] == "-stop": stop = int(sys.argv[i+1]) if start < 0 and stop < 1: print("must specify a start and/or stop step", file=sys.stderr) raise SystemExit data = GetFile(len(sys.argv)-1) if start < 0: start = 0 if stop < 0: stop = data.NumSteps()-1 data.WriteFile(sys.stdout,start,stop) elif sys.argv[1] == "-starts": data = GetFile(2) for i in range(len(data.Steps)): print(i, end=' ') if data.Steps[i].FrontData[0].Flags & CF_HAVE_START: print(data.Steps[i].FrontData[0].StartPoint) else: print("None") elif sys.argv[1] == "-merge": files = [] for i in range(2,len(sys.argv)): files.append(sys.argv[i]) data = GrowthData(files[0]) for i in range(1,len(files)): tmp = GrowthData(files[i]) for j in range(len(tmp.Steps)): data.Steps.append(tmp.Steps[j]) data.WriteFile(sys.stdout) elif sys.argv[1] == "-path": pos = float(sys.argv[2]) data = GetFile(3) for i in range(len(data.Steps)): crds,Ks = FindKAtPos(pos,data.Steps[i].FrontData[0].PointData) print(crds,Ks) elif sys.argv[1] == "-nearest_path": pos = float(sys.argv[2]) data = GetFile(3) sif_calc = SifHistCalc(data) fronts = [(i,0) for i in range(data.NumSteps())] path = sif_calc.ClosestPath(pos,fronts) length = 0 for i in range(len(path)): print(length,path[i].Ks[0]) if i < len(path)-1: length += (path[i+1].CartCoord-path[i].CartCoord).Magnitude() else: print("unrecongnized command: (%s)" % (sys.argv[0],), file=sys.stderr)