nugrid_set.py 文件源码

def write_gce_input_wind_ejection_rate(self,file_name="isotopic_table.txt",final_models=[]):

        '''
            Calculates total mass ejection rate from stellar wind phase.
        '''

        e_kin_array=[]
        mz=[]
        m_final=[]
        mz=[]
        ini_m=[]
        p=-1
        for i in range(len(self.runs_H5_out)):
            p+=1
            sefiles=se(self.runs_H5_out[i])
            sefiles_restart=se(self.runs_H5_restart[i])
            mass=sefiles.get("mini")
            metallicity=sefiles.get("zini") 
            mz1='(M='+str(round(mass,2))+',Z='+str(metallicity)+')'
            endcycle=int(sefiles.se.cycles[-1])
            if len(final_models)>0:
                cycs=self.get_last_cycle([0,final_models[p]+500,500],sefiles,sefiles_restart)
            else:
                cycs=self.get_last_cycle([0,-1,500],sefiles,sefiles_restart)
            endcycle=cycs[-1]
            h_mass_frac=sefiles.get(endcycle,'H-1')
            mass_array=sefiles.get(endcycle,'mass')
            #earlier: 1.e-1 but now as in set1 scripts 0.05
            h_free_core=mass_array[np.where(h_mass_frac<5.e-2)[0][-1]]
            m_final.append(h_free_core)     
            mz.append(mz1)
            ini_m.append(mass)

                f1=open(file_name,'r')
                lines=f1.readlines()
                f1.close()
                i=-1
                line1=''
                while (True):
                        i+=1
                        if i>len(lines)-1:
                                break
                        line=lines[i]
                        line1+=lines[i]
                        for k in range(len(mz)):
                                if mz[k] in lines[i]:
                    ## check ahead for H Mfinal:
                    #for h in range(i,i+10):
                    #   if 'H Lifetime:' in lines[h]:
                    #       lifetime=float(lines[h].split(':')[1])
                    #       break
                    #   if h==i+9:
                    #       print 'ERROR: no lifetime in table found!!'
                    #       return
                    ejection=(ini_m[k]-m_final[k]) #/lifetime
                                        line1+=('H Wind ejection: '+'{:.3E}'.format(ejection)+'\n')
                                        break

                f1=open(file_name,'w')
                f1.write(line1)
                f1.close()