Python constants.value函数代码示例

def e2pz(w1,w2,th):
	"""
	Calculates the momentum scale and the relativistic Compton cross section 
	correction according to P. Holm, PRA 37, 3706 (1988).
	from KH 29.05.96
	input:
	w1 = incident energy  [keV]
	w2 = scattered energy [keV]
	th = scattering angle [deg]
	returns:
	pz = momentum scale   [a.u.]
	cf = cross section correction factor such that:
	J(pz) = cf * d^2(sigma)/d(w2)*d(Omega) [barn/atom/keV/srad]
	"""
	w1  = np.array(w1)    # make sure arrays are used
	w2  = np.array(w2)           
	m   = constants.value('electron mass energy equivalent in MeV')*1e3 #511.003      # Mass in natural units
	th  = math.radians(th) # th/180.0*np.pi  # Angle to radians
	alp = constants.value('fine-structure constant') #1.0/137.036  # Fine structure constant
	r0  = constants.value('classical electron radius') #2.8179e-15   # Electron radius
	q   = np.sqrt(w1**2 + w2**2-2*w1*w2*np.cos(th))                        # Momentum transfer    
	pz  = q/2.0 - (w1-w2) * np.sqrt(1/4.0 + m**2/(2*w1*w2*(1-np.cos(th)))) # In natural units
	E   = np.sqrt(m**2+pz**2)
	A   = ((w1-w2)*E-w1*w2*(1-np.cos(th)))/q
	D   = (w1-w2*np.cos(th))*A/q
	R   = w1*(E-D)
	R2  = R-w1*w2*(1-np.cos(th))
	chi = R/R2 + R2/R + 2.0*m**2 * (1/R-1/R2) + m**4 * (1/R-1/R2)**2
	cf  = 2.0*w1*q*E/(m**2*r0**2*w2*chi)
	cf  = cf*(1e-28*(m*alp))			                             # Cross section now in barns/atom/keV/srad
	pz  = pz/(m*alp)                                                 # pz to atomic units (a.u.)
	return pz, cf

def electrostatic_units(energy):
    """ coefficient to convert to appropriate electrostatic units """
    area = (1 / ((2 * sqrt(2)) ** 2)) * 2 * sqrt(3)
    factor = (1j * ((2 *
              constants.value("Rydberg constant times hc in eV")) ** 5) *
              1e-5 * 2.08e-15 *
            ((constants.value("lattice parameter of silicon") / 1e-10) ** 3))\
              / (area * ((energy) ** 3))
    return factor

    def loadCIF(self,event):
    
        wildcards = 'XRD cifile (*.cif)|*.cif|All files (*.*)|*.*'
        dlg = wx.FileDialog(self, message='Choose CIF',
                           defaultDir=os.getcwd(),
                           wildcard=wildcards, style=wx.FD_OPEN)

        path, read = None, False
        if dlg.ShowModal() == wx.ID_OK:
            read = True
            path = dlg.GetPath().replace('\\', '/')
        dlg.Destroy()
        
        if read:
            cifile = os.path.split(path)[-1]

            try:
                cif = xu.materials.Crystal.fromCIF(path)
            except:
                print('incorrect file format: %s' % os.path.split(path)[-1])
                return

            ## generate hkl list
            hkllist = []
            maxhkl = 8
            for i in range(-maxhkl,maxhkl+1):
                for j in range(-maxhkl,maxhkl+1):
                    for k in range(-maxhkl,maxhkl+1):
                        if i+j+k > 0: # as long as h,k,l all positive, eliminates 0,0,0
                            hkllist.append([i,j,k])
            
            hc = constants.value(u'Planck constant in eV s') * \
                       constants.value(u'speed of light in vacuum') * 1e-3 ## units: keV-m

            if self.wavelength is not None:
                qlist = cif.Q(hkllist)
                Flist = cif.StructureFactorForQ(qlist,(hc/(self.wavelength*(1e-10))*1e3))
            
                Fall = []
                qall = []
                hklall = []
                for i,hkl in enumerate(hkllist):
                    if np.abs(Flist[i]) > 0.01:
                        Fadd = np.abs(Flist[i])
                        qadd = np.linalg.norm(qlist[i])
                        if qadd not in qall and qadd < 6:
                            Fall.extend((0,Fadd,0))
                            qall.extend((qadd,qadd,qadd))
                if np.shape(Fall)[0] > 0:
                    Fall = np.array(Fall)
                    qall = np.array(qall)
                    self.add1Ddata(qall,Fall,name=os.path.split(path)[-1],cif=True)
                else:
                    print('Could not calculate real structure factors.')
            else:
                print('Wavelength/energy must be specified for structure factor calculations.')

 def __init__(self, frep_MHz = None, n = None):
     if frep_MHz is not None:
         self._frep_MHz = frep_MHz
         if frep_MHz > 1.0e6:
             warnings.warn("frep should be specified in MHz; large value given.")
     if n is not None:
         self.set_NPTS(n)            
     # Constants, moved here so that module runs through Sphynx autodoc when
     # scipy is Mocked out.
     self._c_nmps = constants.value('speed of light in vacuum')*1e9/1e12 # c in nm/ps
     self._c_mks  = constants.value('speed of light in vacuum') # m/s        

    def __init__(self):
    
        ## Constructor
        dialog = wx.Dialog.__init__(self, None, title='Define wavelength/energy')#,size=(500, 440))
        ## remember: size=(width,height)
        
        panel = wx.Panel(self)
        
        main = wx.BoxSizer(wx.VERTICAL)
        hmain1 = wx.BoxSizer(wx.HORIZONTAL)
        
        #####
        ## Energy or Wavelength
        hmain1 = wx.BoxSizer(wx.HORIZONTAL)
        self.ch_EorL = wx.Choice(panel,choices=['Energy (keV)','Wavelength (A)'])
        self.entr_EorL = wx.TextCtrl(panel)#, size=(110, -1))
 
        self.ch_EorL.Bind(wx.EVT_CHOICE, self.onEorLSel)
 
        hmain1.Add(self.ch_EorL,   flag=wx.RIGHT,  border=8)
        hmain1.Add(self.entr_EorL, flag=wx.RIGHT,  border=8)
        
        #####
        ## OKAY!
        hmain2 = wx.BoxSizer(wx.HORIZONTAL)
        #hlpBtn = wx.Button(panel, wx.ID_HELP    )
        okBtn  = wx.Button(panel, wx.ID_OK      )
        canBtn = wx.Button(panel, wx.ID_CANCEL  )

        #hmain2.Add(hlpBtn,flag=wx.RIGHT,  border=8)
        hmain2.Add(canBtn, flag=wx.RIGHT, border=8) 
        hmain2.Add(okBtn,  flag=wx.RIGHT,  border=8)

        main.Add(hmain1, flag=wx.ALL, border=10) 
        main.Add(hmain2, flag=wx.ALL|wx.ALIGN_RIGHT, border=10) 

        panel.SetSizer(main)

        self.Show()
        ix,iy = panel.GetBestSize()
        self.SetSize((ix+20, iy+20))
        
        ## set default
        self.hc = constants.value(u'Planck constant in eV s') * \
                       constants.value(u'speed of light in vacuum') * 1e-3 ## units: keV-m
        self.energy = 19.0        
        self.wavelength = self.hc/(self.energy)*1e10 ## units: A
        self.ch_EorL.SetSelection(0)
        self.entr_EorL.SetValue(str(self.energy))

def convert():
    """ converts xyz file from angstroms to bohrs """
    conversion_factor = constants.angstrom / constants.value("Bohr radius")
    xx_coord, yy_coord, zz_coord = genfromtxt(XYZ_CONVERT, unpack="True")
    xyz = column_stack((xx_coord, yy_coord, zz_coord))
    new_xyz = xyz * 1/conversion_factor
    savetxt(XYZ_CONVERT_OUT, new_xyz, fmt=('%3.15f'), delimiter=' ')

 def __init__(self,Ts,qE0T,phis,lamb,length=22.e-3):
     self.Ts       = Ts
     self.Tsf      = Ts     # final kin. energy
     self.qE0T     = qE0T
     self.phis     = phis   # rad
     self.lamb     = lamb
     self.length   = length
     self.nbslices = 10     # nbof slices
     dl            = length/self.nbslices
     m0c2          = C.value('proton mass energy equivalent in MeV')  # MeV
     g             = 1.+Ts/m0c2   # Lorentz gamma
     bg            = sqrt(g**2-1.)
     beta          = bg/g         # Lorentz beta
     # phases,energies & maps of slices
     phases    = []
     dtks      = []
     self.maps = []
     for i in range(self.nbslices+1):
         phase = self.phis+2*pi/(beta*lamb)*(i*dl-self.length/2.)
         phases.append(phase)
         dtk = qE0T*dl*cos(phase)
         dtks.append(dtk)
     for i in range(self.nbslices):
         dgdp,f,pot,ham = Hamiltonian(self.Tsf,qE0T,phases[i],lamb)
         self.maps.append(Order3Map(+dl,dgdp,f))
         self.Tsf += dtks[i]
     self.Tsf -= dtks[i]     # one too much

def Hamiltonian(Ts,qE0T,phis,lamb):       
    """ 
        kanonische Koordinaten: p (aka w), x (aka phi)
        unabhaegige Variable: t (aka s)
        Hamiltonian H(p,x,t) = g(p) + V(x,t)
        g = A*p**2/2
        V = B*(sin(x)-x*cos(phis))
        f = -dV/dx = -B*(cos(x)-cos(phis))
        dgdp = dH/dp = A*p
    """
    # closure
    m0c2 = C.value('proton mass energy equivalent in MeV')  # MeV
    g    = 1+Ts/m0c2
    bgs  = sqrt(g**2-1.)
    A    = 2.*pi/bgs**3/lamb
    B    = qE0T/m0c2
    # potential
    def V(x,B=B,phis=phis):
        return B*(sin(x)-x*cos(phis))
    # hamiltonian
    def H(p,x,B=B,phis=phis):
        return A*p**2+V(x)
    # dH/dp
    def dgdp(p,A=A):
        return A*p
    # -dV/dx
    def f(x,t,B=B,phis=phis):
        return -B*(cos(x)-cos(phis))

    return dgdp,f,V,H

def freq_units(units):
    """
    Returns conversion factor from THz to the requred units and the label in the form of a namedtuple
    Accepted values: thz, ev, mev, ha, cm-1, cm^-1
    """

    d = {"thz": FreqUnits(1, "THz"),
         "ev": FreqUnits(const.value("hertz-electron volt relationship") * const.tera, "eV"),
         "mev": FreqUnits(const.value("hertz-electron volt relationship") * const.tera / const.milli, "meV"),
         "ha": FreqUnits(const.value("hertz-hartree relationship") * const.tera, "Ha"),
         "cm-1": FreqUnits(const.value("hertz-inverse meter relationship") * const.tera * const.centi, "cm^{-1}"),
         'cm^-1': FreqUnits(const.value("hertz-inverse meter relationship") * const.tera * const.centi, "cm^{-1}")
         }
    try:
        return d[units.lower().strip()]
    except KeyError:
        raise KeyError('Value for units `{}` unknown\nPossible values are:\n {}'.format(units, list(d.keys())))

 def rotational_constants(self, units='ghz'):
     """Compute the rotational constants in 1/cm or GHz."""
     choices = ('invcm', 'ghz')
     units = units.lower()
     if units not in choices:
         raise ValueError("Invalid units, pick one of {}".format(choices))
     principal_moments = self.principal_moments_of_inertia()[0]
     _check_scipy(_found_scipy)
     bohr2ang = spc.value('atomic unit of length') / spc.angstrom
     xfamu = 1 / spc.value('electron mass in u')
     xthz = spc.value('hartree-hertz relationship')
     rotghz = xthz * (bohr2ang ** 2) / (2 * xfamu * spc.giga)
     if units == 'ghz':
         conv = rotghz
     if units == 'invcm':
         ghz2invcm = spc.giga * spc.centi / spc.c
         conv = rotghz * ghz2invcm
     return conv / principal_moments

def water_term_pol(B0, T):
    gamma = constants.value('proton gyromag. ratio')
    hbar = constants.hbar
    k = constants.k

    larmor_w = gamma * B0
    pol = hbar * larmor_w
    pol /= (2 * k * T)
    return pol

    def _read_funcfl_file(self):
        rphi = np.sqrt(27.2*0.529)
        eV2J = spc.value('electron volt')

        print('Reading EAM potential in single-element (funcfl) format')
        f = open(self.filename, 'r')

        line = f.readline()  # comment line
        self.atomic_number = int(f.readline().split()[0])

        line = f.readline().split()
        self.nrho    =   int(line[0]) 
        self.drho    = float(line[1])
        self.nr      =   int(line[2])
        self.dr      = float(line[3])
        self.cutoff  = float(line[4])

        self.fembed  = np.zeros((self.nrho,2))
        self.effchg  = np.zeros((self.nr,2))
        self.fdens   = np.zeros((self.nr,2))

        # read embedding functional
        line    = f.readline().split()
        n_col   = len(line)
        n_lines = int(self.nrho/self.ncol)-1

        i_rows = 0
        while i_rows < self.nrho:
            for val in line:
                self.fembed[i_rows,0] = i_rows*self.drho
                self.fembed[i_rows,1] = float(val)
                i_rows += 1
                if i_rows >= self.nrho: 
                    break
            line=f.readline().split()
        
        # read charge function
        i_rows = 0
        while i_rows < self.nr:
            for val in line:
                self.effchg[i_rows,0]=i_rows*self.dr
                self.effchg[i_rows,1]=float(val)*rphi/(1.e-30+self.effchg[self.i_rows,0])
                i_rows += 1
                if i_rows >= self.nr: 
                    break
            line=f.readline().split()

        # read density function
        i_rows = 0
        while i_rows < self.nr:
            for val in line:
                self.fdens[i_rows,0]=i_rows*self.dr
                self.fdens[i_rows,1]=float(val)
                i_rows += 1
                if i_rows >= self.nr:
                    break
            line=f.readline().split()

    def __init__(self, N=DEFAULT_N):
        """
        Parameters
        ----------
        N : int
            the number of grid points, it should be a power of 2

        Examples
        --------
        >>> from generic_potential import GenericPotential
        >>> device = GenericPotential(1024)
        """

        # default grid size
        self.N = N
        self.hN = int(self.N/2) # half grid
        self.points_before = int(self.N/2)
        self.points_after = int(self.N/2)

        # useful atomic unities
        self.au_l    = cte.value('atomic unit of length')
        self.au_t    = cte.value('atomic unit of time')
        self.au_e    = cte.value('atomic unit of energy')
        self.au_v    = cte.value('atomic unit of electric potential')
        self.hbar_au = 1.0
        self.me_au   = 1.0

        # other useful constants
        self.ev = cte.value('electron volt')
        self.c  = cte.value('speed of light in vacuum') # m/s
        self.me = cte.value('electron mass')
        self.q  = cte.value('elementary charge')

        # relations of interest
        self.au2ev  = self.au_e / self.ev
        self.au2ang = self.au_l / 1e-10

        # specific for default quantum harmonic oscillator
        self.l     = 0.0000081 # m
        self.f     = self.c / self.l # Hz
        self.w     = 2.0 * np.pi * self.f
        self.z_m   = np.linspace(-5e-9,5e-9, self.N)
        self.v_j   = 0.5 * self.me * self.z_m**2 * self.w**2
        self.z_nm  = self.z_m / 1e-9
        self.v_ev  = self.v_j / self.ev
        self.m_eff = np.ones(self.z_nm.size)

        # set default time increment
        self._set_dt()

        # adjust grid
        self.normalize_device()

    def onLoadCIF(self, event=None):

        wildcards = 'CIF file (*.cif)|*.cif|All files (*.*)|*.*'
        dlg = wx.FileDialog(self, message='Choose CIF file',
                           defaultDir=os.getcwd(),
                           wildcard=wildcards, style=wx.FD_OPEN)

        path, read = None, False
        if dlg.ShowModal() == wx.ID_OK:
            read = True
            path = dlg.GetPath().replace('\\', '/')
        dlg.Destroy()

        if read:
            cry_strc = struc_from_cif(path)

            if cry_strc:
                hc = constants.value(u'Planck constant in eV s') * \
                         constants.value(u'speed of light in vacuum') * 1e-3 ## units: keV-m
                energy = hc/(self.xrd.wavelength) ## units: keV
                q,F = calc_all_F(cry_strc,energy,maxhkl=10,qmax=5)

                #self.plot1Dxrd(xrddata, show_xrd2=True)
                print('Values are calculated; plotting not yet implemented.')

 def principal_moments_of_inertia(self, units='amu_bohr_2'):
     """Return the principal moments of inertia in 3 kinds of units:
     1. [amu][bohr]^2
     2. [amu][angstrom]^2
     3. [g][cm]^2
     and the principal axes.
     """
     choices = ('amu_bohr_2', 'amu_angstrom_2', 'g_cm_2')
     units = units.lower()
     if units not in choices:
         raise ValueError("Invalid units, pick one of {}".format(choices))
     moi_tensor = self.moment_of_inertia_tensor()
     principal_moments, principal_axes = np.linalg.eigh(moi_tensor)
     if units == 'amu_bohr_2':
         conv = 1
     if units == 'amu_angstrom_2':
         _check_scipy(_found_scipy)
         bohr2ang = spc.value('atomic unit of length') / spc.angstrom
         conv = bohr2ang ** 2
     if units == 'g_cm_2':
         _check_scipy(_found_scipy)
         amu2g = spc.value('unified atomic mass unit') * spc.kilo
         conv = amu2g * (spc.value('atomic unit of length') * spc.centi) ** 2
     return conv * principal_moments, principal_axes