\(\renewcommand\AA{\text{Å}}\)
12. GSASIImath: computation module
12.1. Summary/Contents
Least-squares minimization and misc. computation routines.
12.2. GSASIImath Classes and routines
Routines defined in GSASIImath
follow.
- GSASIImath.ApplyModeDisp(data)[source]
Applies ISODISTORT mode displacements to atom lists. This changes the contents of the Draw Atoms positions and the Atoms positions.
- Parameters:
data (dict) – the contents of the Phase data tree item for a particular phase
- GSASIImath.ApplyModulation(data, tau)[source]
Applies modulation to drawing atom positions & Uijs for given tau
- GSASIImath.ApplySeqData(data, seqData, PF2=False)[source]
Applies result from seq. refinement to drawing atom positions & Uijs
- Parameters:
data (dict) – GSAS-II phase data structure
seqData (dict) – GSAS-II sequential refinement results structure
PF2 (bool) – if True then seqData is from a sequential run of PDFfit2
- Returns:
list drawAtoms: revised Draw Atoms list
- GSASIImath.AtomTLS2UIJ(atomData, atPtrs, Amat, rbObj)[source]
default doc string
- Parameters:
name (type) – description
- Returns:
type name: description
- GSASIImath.AtomsCollect(data, Ind, Sel)[source]
Finds the symmetry set of atoms for those selected. Selects the one closest to the selected part of the unit cell for the selected atoms
- Parameters:
data – the phase data structure
Ind (list) – list of selected atom indices
Sel (int) – an index with the selected plane or location in the unit cell to find atoms closest to
- Returns:
the list of unique atoms where selected atoms may have been replaced. Anisotropic Uij’s are transformed
- GSASIImath.BessIn(nmax, x)[source]
compute modified Bessel function I(n,x) from scipy routines & recurrance relation returns sequence of I(n,x) for n in range [-nmax…0…nmax]
- Parameters:
nmax (integer) – maximul order for In(x)
x (float) – argument for In(x)
- Returns numpy array:
[I(-nmax,x)…I(0,x)…I(nmax,x)]
- GSASIImath.BessJn(nmax, x)[source]
compute Bessel function J(n,x) from scipy routine & recurrance relation returns sequence of J(n,x) for n in range [-nmax…0…nmax]
- Parameters:
nmax (integer) – maximul order for Jn(x)
x (float) – argument for Jn(x)
- Returns numpy array:
[J(-nmax,x)…J(0,x)…J(nmax,x)]
- GSASIImath.CalcAngle(angle_dict, angle_atoms, parmDict)[source]
Used in class:GSASIIobj.ExpressionCalcObj to compute bond angles when defined in an expression.
- GSASIImath.CalcAngleDeriv(angle_dict, angle_atoms, parmDict)[source]
Used to compute s.u. values on angles tracked in the sequential results table
- GSASIImath.CalcDist(distance_dict, distance_atoms, parmDict)[source]
Used in class:GSASIIobj.ExpressionCalcObj to compute bond distances when defined in an expression.
- GSASIImath.CalcDistDeriv(distance_dict, distance_atoms, parmDict)[source]
Used to compute s.u. values on distances tracked in the sequential results table
- GSASIImath.CalcIsoCoords(Phase, parmDict, covdata={})[source]
Compute the coordinate positions from ISODISTORT displacement mode values Uncertainties are computed if covdata is supplied.
- Parameters:
Phase (dict) – contents of tree entry for selected phase
parmDict (dict) – a dict with values for the modes; note that in the parmDict from refinements the mode values are not normalized, but this assumes they are.
Phase – full covariance information from tree
- Returns:
modeDict,posDict where modeDict contains pairs of mode values and mode s.u. values; posDict contains pairs of displacement values and their s.u. values.
- GSASIImath.CalcIsoDisp(Phase, parmDict={}, covdata={})[source]
Compute the ISODISTORT displacement variable values from the atomic coordinates, applying the p::dA?:n displacements if parmDict is supplied. Uncertainties are computed if covdata is supplied.
- GSASIImath.ChargeFlip(data, reflDict, pgbar)[source]
default doc string
- Parameters:
name (type) – description
- Returns:
type name: description
- GSASIImath.Den2Vol(Elements, density)[source]
converts density to molecular volume
- Parameters:
Elements (dict) – elements in molecular formula; each must contain Num: number of atoms in formula Mass: at. wt.
density (float) – material density in gm/cm^3
- Returns:
float volume: molecular volume in A^3
- GSASIImath.DrawAtomsReplaceByID(data, loc, atom, ID)[source]
Replace all atoms in drawing array with an ID matching the specified value
- GSASIImath.El2EstVol(Elements)[source]
Estimate volume from molecular formula; assumes atom volume = 10A^3
- Parameters:
Elements (dict) – elements in molecular formula; each must contain Num: number of atoms in formula
- Returns:
float volume: estimate of molecular volume in A^3
- GSASIImath.El2Mass(Elements)[source]
compute molecular weight from Elements
- Parameters:
Elements (dict) – elements in molecular formula; each must contain Num: number of atoms in formula Mass: at. wt.
- Returns:
float mass: molecular weight.
- GSASIImath.FillAtomLookUp(atomData, indx)[source]
create a dictionary of atom indexes with atom IDs as keys
- Parameters:
atomData (list) – Atom table to be used
indx (int) – pointer to position of atom id in atom record (typically cia+8)
- Returns:
dict atomLookUp: dictionary of atom indexes with atom IDs as keys
- GSASIImath.FindAllNeighbors(phase, FrstName, AtNames, notName='', Orig=None, Short=False, searchType='Bond')[source]
Find neighboring atoms Uses Bond search criteria unless searchType is set to non-default
- GSASIImath.FindAtomIndexByIDs(atomData, loc, IDs, Draw=True)[source]
finds the set of atom array indices for a list of atom IDs. Will search either the Atom table or the drawAtom table.
- Parameters:
atomData (list) – Atom or drawAtom table containting coordinates, etc.
loc (int) – location of atom id in atomData record
IDs (list) – atom IDs to be found
Draw (bool) – True if drawAtom table to be searched; False if Atom table is searched
- Returns:
list indx: atom (or drawAtom) indices
- GSASIImath.Fourier4DMap(data, reflDict)[source]
default doc string
- Parameters:
name (type) – description
- Returns:
type name: description
- GSASIImath.FourierMap(data, reflDict)[source]
default doc string
- Parameters:
name (type) – description
- Returns:
type name: description
- GSASIImath.GetAngleSig(Oatoms, Atoms, Amat, SGData, covData={})[source]
Not used
- Parameters:
name (type) – description
- Returns:
type name: description
- GSASIImath.GetAtomCoordsByID(pId, parmDict, AtLookup, indx)[source]
default doc string
- Parameters:
name (type) – description
- Returns:
type name: description
- GSASIImath.GetAtomFracByID(pId, parmDict, AtLookup, indx)[source]
default doc string
- Parameters:
name (type) – description
- Returns:
type name: description
- GSASIImath.GetAtomItemsById(atomData, atomLookUp, IdList, itemLoc, numItems=1)[source]
gets atom parameters for atoms using atom IDs
- Parameters:
atomData (list) – Atom table to be used
atomLookUp (dict) – dictionary of atom indexes with atom IDs as keys
IdList (list) – atom IDs to be found
itemLoc (int) – pointer to desired 1st item in an atom table entry
numItems (int) – number of items to be retrieved
- Returns:
type name: description
- GSASIImath.GetAtomMomsByID(pId, parmDict, AtLookup, indx)[source]
default doc string
- Parameters:
name (type) – description
- Returns:
type name: description
- GSASIImath.GetAtomsById(atomData, atomLookUp, IdList)[source]
gets a list of atoms from Atom table that match a set of atom IDs
- Parameters:
atomData (list) – Atom table to be used
atomLookUp (dict) – dictionary of atom indexes with atom IDs as keys
IdList (list) – atom IDs to be found
- Returns:
list atoms: list of atoms found
- GSASIImath.GetDATSig(Oatoms, Atoms, Amat, SGData, covData={})[source]
default doc string
- Parameters:
name (type) – description
- Returns:
type name: description
- GSASIImath.GetDistSig(Oatoms, Atoms, Amat, SGData, covData={})[source]
not used
- Parameters:
name (type) – description
- Returns:
type name: description
- GSASIImath.GetMag(mag, Cell)[source]
Compute magnetic moment magnitude. :param list mag: atom magnetic moment parms (must be magnetic!) :param list Cell: lattice parameters
- Returns:
moment magnitude as float
- GSASIImath.GetMagDerv(mag, Cell)[source]
Compute magnetic moment derivatives numerically :param list mag: atom magnetic moment parms (must be magnetic!) :param list Cell: lattice parameters
- Returns:
moment derivatives as floats
- GSASIImath.GetSHCoeff(pId, parmDict, SHkeys)[source]
default doc string
- Parameters:
name (type) – description
- Returns:
type name: description
- GSASIImath.GetTorsionSig(Oatoms, Atoms, Amat, SGData, covData={})[source]
default doc string
- Parameters:
name (type) – description
- Returns:
type name: description
- GSASIImath.GetXYZDist(xyz, XYZ, Amat)[source]
- gets distance from position xyz to all XYZ, xyz & XYZ are np.array
and are in crystal coordinates; Amat is crystal to Cart matrix
- Parameters:
name (type) – description
- Returns:
type name: description
- GSASIImath.HessianLSQ(func, x0, Hess, args=(), ftol=1.49012e-08, xtol=1e-06, maxcyc=0, lamda=-3, Print=False, refPlotUpdate=None)[source]
Minimize the sum of squares of a function (\(f\)) evaluated on a series of values (y): \(\sum_{y=0}^{N_{obs}} f(y,{args})\) where \(x = arg min(\sum_{y=0}^{N_{obs}} (func(y)^2,axis=0))\)
- Parameters:
func (function) – callable method or function should take at least one (possibly length N vector) argument and returns M floating point numbers.
x0 (np.ndarray) – The starting estimate for the minimization of length N
Hess (function) – callable method or function A required function or method to compute the weighted vector and Hessian for func. It must be a symmetric NxN array
args (tuple) – Any extra arguments to func are placed in this tuple.
ftol (float) – Relative error desired in the sum of squares.
xtol (float) – Relative tolerance of zeros in the SVD solution in nl.pinv.
maxcyc (int) – The maximum number of cycles of refinement to execute, if -1 refine until other limits are met (ftol, xtol)
lamda (int) – initial Marquardt lambda=10**lamda
Print (bool) – True for printing results (residuals & times) by cycle
- Returns:
(x,cov_x,infodict) where
x : ndarray The solution (or the result of the last iteration for an unsuccessful call).
cov_x : ndarray Uses the fjac and ipvt optional outputs to construct an estimate of the jacobian around the solution. This matrix must be multiplied by the residual standard deviation to get the covariance of the parameter estimates – see curve_fit.
– or
None
if a singular matrix encountered (indicates very flat curvature in direction), or some other error is encountered.– or an empty array if a refinement was not performed because no valid variables were refined.
infodict : dict, a dictionary of optional outputs with the keys:
’fvec’ : the function evaluated at the output
’num cyc’:
’nfev’: number of objective function evaluation calls
’lamMax’:
’psing’: list of variable variables that have been removed from the refinement
’SVD0’: -1 for singlar matrix, -2 for objective function exception, Nzeroes = # of SVD 0’s
’Hcorr’: list entries (i,j,c) where i & j are of highly correlated variables & c is correlation coeff.
- GSASIImath.HessianSVD(func, x0, Hess, args=(), ftol=1.49012e-08, xtol=1e-06, maxcyc=0, lamda=-3, Print=False, refPlotUpdate=None)[source]
Minimize the sum of squares of a function (\(f\)) evaluated on a series of values (y): \(\sum_{y=0}^{N_{obs}} f(y,{args})\) where \(x = arg min(\sum_{y=0}^{N_{obs}} (func(y)^2,axis=0))\)
- Parameters:
func (function) – callable method or function should take at least one (possibly length N vector) argument and returns M floating point numbers.
x0 (np.ndarray) – The starting estimate for the minimization of length N
Hess (function) – callable method or function A required function or method to compute the weighted vector and Hessian for func. It must be a symmetric NxN array
args (tuple) – Any extra arguments to func are placed in this tuple.
ftol (float) – Relative error desired in the sum of squares.
xtol (float) – Relative tolerance of zeros in the SVD solution in nl.pinv.
maxcyc (int) – The maximum number of cycles of refinement to execute, if -1 refine until other limits are met (ftol, xtol)
Print (bool) – True for printing results (residuals & times) by cycle
- Returns:
(x,cov_x,infodict) where
x : ndarray The solution (or the result of the last iteration for an unsuccessful call).
cov_x : ndarray Uses the fjac and ipvt optional outputs to construct an estimate of the jacobian around the solution.
None
if a singular matrix encountered (indicates very flat curvature in some direction). This matrix must be multiplied by the residual standard deviation to get the covariance of the parameter estimates – see curve_fit.infodict : dict a dictionary of optional outputs with the keys:
’fvec’ : the function evaluated at the output
’num cyc’:
’nfev’:
’lamMax’:0.
’psing’:
’SVD0’:
- GSASIImath.MagMod(glTau, xyz, modQ, MSSdata, SGData, SSGData)[source]
this needs to make magnetic moment modulations & magnitudes as fxn of gTau points; NB: this allows only 1 mag. wave fxn.
- GSASIImath.Modulation(H, HP, nWaves, Fmod, Xmod, Umod, glTau, glWt)[source]
H: array nRefBlk x ops X hklt HP: array nRefBlk x ops X hklt proj to hkl nWaves: list number of waves for frac, pos, uij & mag Fmod: array 2 x atoms x waves (sin,cos terms) Xmod: array atoms X 3 X ngl Umod: array atoms x 3x3 x ngl glTau,glWt: arrays Gauss-Lorentzian pos & wts
- GSASIImath.ModulationDerv(H, HP, Hij, nWaves, waveShapes, Fmod, Xmod, UmodAB, SCtauF, SCtauX, SCtauU, glTau, glWt)[source]
Compute Fourier modulation derivatives H: array ops X hklt proj to hkl HP: array ops X hklt proj to hkl Hij: array 2pi^2[a*^2h^2 b*^2k^2 c*^2l^2 a*b*hk a*c*hl b*c*kl] of projected hklm to hkl space
- GSASIImath.ModulationTw(H, HP, nWaves, Fmod, Xmod, Umod, glTau, glWt)[source]
H: array nRefBlk x tw x ops X hklt HP: array nRefBlk x tw x ops X hklt proj to hkl Fmod: array 2 x atoms x waves (sin,cos terms) Xmod: array atoms X ngl X 3 Umod: array atoms x ngl x 3x3 glTau,glWt: arrays Gauss-Lorentzian pos & wts
- GSASIImath.NCScattDen(Elements, vol, wave=0.0)[source]
Estimate neutron scattering density from molecular formula & volume; ignores valence, but includes anomalous effects
- Parameters:
Elements (dict) – elements in molecular formula; each element must contain Num: number of atoms in formula Z: atomic number
vol (float) – molecular volume in A^3
wave (float) – optional wavelength in A
- Returns:
float rho: scattering density in 10^10cm^-2; if wave > 0 the includes f’ contribution
- Returns:
float mu: if wave>0 absorption coeff in cm^-1 ; otherwise 0
- Returns:
float fpp: if wave>0 f” in 10^10cm^-2; otherwise 0
- GSASIImath.OmitMap(data, reflDict, pgbar=None)[source]
default doc string
- Parameters:
name (type) – description
- Returns:
type name: description
- GSASIImath.PeaksEquiv(data, Ind)[source]
Find the equivalent map peaks for those selected. Works on the contents of data[‘Map Peaks’].
- Parameters:
data – the phase data structure
Ind (list) – list of selected peak indices
- Returns:
augmented list of peaks including those related by symmetry to the ones in Ind
- GSASIImath.PeaksUnique(data, Ind, Sel, dlg)[source]
Finds the symmetry unique set of peaks from those selected. Selects the one closest to the center of the unit cell. Works on the contents of data[‘Map Peaks’]. Called from OnPeaksUnique in GSASIIphsGUI.py,
- Parameters:
data – the phase data structure
Ind (list) – list of selected peak indices
Sel (int) – selected column to find peaks closest to
dlg (wx object) – progress bar dialog box
- Returns:
the list of symmetry unique peaks from among those given in Ind
- GSASIImath.Q2AV(Q)[source]
convert quaternion to angle (radians 0-2pi) & normalized vector q=r+ai+bj+ck
- GSASIImath.Q2AVdeg(Q)[source]
convert quaternion to angle (degrees 0-360) & normalized vector q=r+ai+bj+ck
- GSASIImath.RotateRBXYZ(Bmat, Cart, oriQ, symAxis=None)[source]
rotate & transform cartesian coordinates to crystallographic ones no translation applied. To be used for numerical derivatives
- Parameters:
Bmat (array) – Orthogonalization matrix, see
GSASIIlattice.cell2AB()
Cart (array) – 2D array of coordinates
Q (array) – quaternion as an np.array
symAxis (tuple) – if not None (default), specifies the symmetry axis of the rigid body, which will be aligned to the quaternion vector.
- Returns:
2D array of fractional coordinates, without translation to origin
- GSASIImath.SSChargeFlip(data, reflDict, pgbar)[source]
default doc string
- Parameters:
name (type) – description
- Returns:
type name: description
- GSASIImath.SearchMap(generalData, drawingData, Neg=False)[source]
Does a search of a density map for peaks meeting the criterion of peak height is greater than mapData[‘cutOff’]/100 of mapData[‘rhoMax’] where mapData is data[‘General’][‘mapData’]; the map is also in mapData.
- Parameters:
generalData – the phase data structure; includes the map
drawingData – the drawing data structure
Neg – if True then search for negative peaks (i.e. H-atoms & neutron data)
- Returns:
(peaks,mags,dzeros) where
peaks : ndarray x,y,z positions of the peaks found in the map
mags : ndarray the magnitudes of the peaks
dzeros : ndarray the distance of the peaks from the unit cell origin
dcent : ndarray the distance of the peaks from the unit cell center
- GSASIImath.SetDefaultDData(dType, histoName, NShkl=0, NDij=0)[source]
Sets default values for various histogram parameters param: str dType: 3 letter histogram type, e.g. ‘PNT’ param: str histoName: histogram name as it aoears in tree param: NShkl int: number of generalized mustrain coefficients - depends on symmetry param: NDij int: number of hydrostatic strain coefficients - depends on symmetry
returns dict: default data for histogram - found in data tab for phase/histogram
- GSASIImath.SetMolCent(model, RBData)[source]
default doc string
- Parameters:
name (type) – description
- Returns:
type name: description
- GSASIImath.TLS2Uij(xyz, g, Amat, rbObj)[source]
default doc string
- Parameters:
name (type) – description
- Returns:
type name: description
- GSASIImath.UpdateHKLFvals(histoName, phaseData, reflData)[source]
Update the flag field and the d-space field in HKLF reflection table.
Data tree contents are passed into this routine as data objects (strings, dicts,…) rather than data tree references so that this can called both from GUI routines as well as used in scripting.
This gets called from
GSASIIphsGUI.CheckAddHKLF()
(which gets called in OnHklfAdd, insideGSASIIphsGUI.UpdatePhaseData()
andGSASIIddataGUI.MakeHistPhaseWin()
). Also, inGSASIIscriptable.G2Project.link_histogram_phase()
and in routines OnImportPhase or OnImportSfact inside func:GSASIIdataGUI.GSASIImain.
- GSASIImath.UpdateMCSAxyz(Bmat, MCSA)[source]
default doc string
- Parameters:
name (type) – description
- Returns:
type name: description
- GSASIImath.UpdateRBUIJ(Bmat, Cart, RBObj)[source]
default doc string
- Parameters:
name (type) – description
- Returns:
type name: description
- GSASIImath.UpdateRBXYZ(Bmat, RBObj, RBData, RBType)[source]
returns crystal coordinates for atoms described by RBObj. Note that RBObj[‘symAxis’], if present, determines the symmetry axis of the rigid body, which will be aligned to the quaternion direction.
- Parameters:
Bmat (np.array) – see
GSASIIlattice.cell2AB()
rbObj (dict) – rigid body selection/orientation information
RBData (dict) – rigid body tree data structure
RBType (str) – rigid body type, ‘Vector’ or ‘Residue’
- Returns:
coordinates for rigid body as XYZ,Cart where XYZ is the location in crystal coordinates and Cart is in cartesian
- GSASIImath.ValEsd(value, esd=0, nTZ=False)[source]
Format a floating point number with a given level of precision or with in crystallographic format with a “esd”, as value(esd). If esd is negative the number is formatted with the level of significant figures appropriate if abs(esd) were the esd, but the esd is not included. if the esd is zero, approximately 6 significant figures are printed. nTZ=True causes “extra” zeros to be removed after the decimal place. for example:
“1.235(3)” for value=1.2346 & esd=0.003
“1.235(3)e4” for value=12346. & esd=30
“1.235(3)e6” for value=0.12346e7 & esd=3000
“1.235” for value=1.2346 & esd=-0.003
“1.240” for value=1.2395 & esd=-0.003
“1.24” for value=1.2395 & esd=-0.003 with nTZ=True
“1.23460” for value=1.2346 & esd=0.0
- Parameters:
value (float) – number to be formatted
esd (float) – uncertainty or if esd < 0, specifies level of precision to be shown e.g. esd=-0.01 gives 2 places beyond decimal
nTZ (bool) – True to remove trailing zeros (default is False)
- Returns:
value(esd) or value as a string
- GSASIImath.Vol2Den(Elements, volume)[source]
converts volume to density
- Parameters:
Elements (dict) – elements in molecular formula; each must contain Num: number of atoms in formula Mass: at. wt.
volume (float) – molecular volume in A^3
- Returns:
float density: material density in gm/cm^3
- GSASIImath.XScattDen(Elements, vol, wave=0.0)[source]
Estimate X-ray scattering density from molecular formula & volume; ignores valence, but includes anomalous effects
- Parameters:
Elements (dict) – elements in molecular formula; each element must contain Num: number of atoms in formula Z: atomic number
vol (float) – molecular volume in A^3
wave (float) – optional wavelength in A
- Returns:
float rho: scattering density in 10^10cm^-2; if wave > 0 the includes f’ contribution
- Returns:
float mu: if wave>0 absorption coeff in cm^-1 ; otherwise 0
- Returns:
float fpp: if wave>0 f” in 10^10cm^-2; otherwise 0
- GSASIImath.adjHKLmax(SGData, Hmax)[source]
default doc string
- Parameters:
name (type) – description
- Returns:
type name: description
- GSASIImath.anneal(func, x0, args=(), schedule='fast', T0=None, Tf=1e-12, maxeval=None, maxaccept=None, maxiter=400, feps=1e-06, quench=1.0, c=1.0, lower=-100, upper=100, dwell=50, slope=0.9, ranStart=False, ranRange=0.1, autoRan=False, dlg=None)[source]
Minimize a function using simulated annealing.
Schedule is a schedule class implementing the annealing schedule. Available ones are ‘fast’, ‘cauchy’, ‘boltzmann’
- Parameters:
func (callable) – f(x, *args) Function to be optimized.
x0 (ndarray) – Initial guess.
args (tuple) – Extra parameters to func.
schedule (base_schedule) – Annealing schedule to use (a class).
T0 (float) – Initial Temperature (estimated as 1.2 times the largest cost-function deviation over random points in the range).
Tf (float) – Final goal temperature.
maxeval (int) – Maximum function evaluations.
maxaccept (int) – Maximum changes to accept.
maxiter (int) – Maximum cooling iterations.
feps (float) – Stopping relative error tolerance for the function value in last four coolings.
quench,c (float) – Parameters to alter fast_sa schedule.
lower,upper (float/ndarray) – Lower and upper bounds on x.
dwell (int) – The number of times to search the space at each temperature.
slope (float) – Parameter for log schedule
ranStart=False (bool) – True for set 10% of ranges about x
- Returns:
(xmin, Jmin, T, feval, iters, accept, retval) where
xmin (ndarray): Point giving smallest value found.
Jmin (float): Minimum value of function found.
T (float): Final temperature.
feval (int): Number of function evaluations.
iters (int): Number of cooling iterations.
accept (int): Number of tests accepted.
retval (int): Flag indicating stopping condition:
0: Points no longer changing
1: Cooled to final temperature
2: Maximum function evaluations
3: Maximum cooling iterations reached
4: Maximum accepted query locations reached
5: Final point not the minimum amongst encountered points
Notes: Simulated annealing is a random algorithm which uses no derivative information from the function being optimized. In practice it has been more useful in discrete optimization than continuous optimization, as there are usually better algorithms for continuous optimization problems.
Some experimentation by trying the difference temperature schedules and altering their parameters is likely required to obtain good performance.
The randomness in the algorithm comes from random sampling in numpy. To obtain the same results you can call numpy.random.seed with the same seed immediately before calling scipy.optimize.anneal.
We give a brief description of how the three temperature schedules generate new points and vary their temperature. Temperatures are only updated with iterations in the outer loop. The inner loop is over range(dwell), and new points are generated for every iteration in the inner loop. (Though whether the proposed new points are accepted is probabilistic.)
For readability, let d denote the dimension of the inputs to func. Also, let x_old denote the previous state, and k denote the iteration number of the outer loop. All other variables not defined below are input variables to scipy.optimize.anneal itself.
In the ‘fast’ schedule the updates are
u ~ Uniform(0, 1, size=d) y = sgn(u - 0.5) * T * ((1+ 1/T)**abs(2u-1) -1.0) xc = y * (upper - lower) x_new = x_old + xc T_new = T0 * exp(-c * k**quench)
- GSASIImath.calcRamaEnergy(phi, psi, Coeff=[])[source]
Computes pseudo potential energy from a pair of torsion angles and a numerical description of the potential energy surface. Used to create penalty function in LS refinement: \(Eval(\phi,\psi) = C[0]*exp(-V/1000)\)
where \(V = -C[3] * (\phi-C[1])^2 - C[4]*(\psi-C[2])^2 - 2*(\phi-C[1])*(\psi-C[2])\)
- Parameters:
phi (float) – first torsion angle (\(\phi\))
psi (float) – second torsion angle (\(\psi\))
Coeff (list) – pseudo potential coefficients
- Returns:
list (sum,Eval): pseudo-potential difference from minimum & value; sum is used for penalty function.
- GSASIImath.calcTorsionEnergy(TOR, Coeff=[])[source]
default doc string
- Parameters:
name (type) – description
- Returns:
type name: description
- GSASIImath.dropTerms(bad, hessian, indices, *vectors)[source]
Remove the ‘bad’ terms from the Hessian and vector
- Parameters:
bad (tuple) – a list of variable (row/column) numbers that should be removed from the hessian and vector. Example: (0,3) removes the 1st and 4th column/row
hessian (np.array) – a square matrix of length n x n
indices (np.array) – the indices of the least-squares vector of length n referenced to the initial variable list; as this routine is called multiple times, more terms may be removed from this list
additional-args – various least-squares model values, length n
- Returns:
hessian, indices, vector0, vector1,… where the lengths are now n’ x n’ and n’, with n’ = n - len(bad)
- GSASIImath.findOffset(SGData, A, Fhkl)[source]
default doc string
- Parameters:
name (type) – description
- Returns:
type name: description
- GSASIImath.findSSOffset(SGData, SSGData, A, Fhklm)[source]
default doc string
- Parameters:
name (type) – description
- Returns:
type name: description
- GSASIImath.fmtPhaseContents(compdict)[source]
Format results from
phaseContents()
- GSASIImath.getAngSig(VA, VB, Amat, SGData, covData={})[source]
Compute an interatomic angle and its uncertainty from two vectors each between an orgin atom and either of a pair of nearby atoms.
- Parameters:
VA (np.array) – an interatomic vector as a structure
VB (np.array) – an interatomic vector also as a structure
Amat (np.array) – unit cell parameters as an A vector
SGData (dict) – symmetry information
covData (dict) – covariance information including the covariance matrix and the list of varied parameters. If not supplied, the s.u. values are returned as zeros.
- Returns:
angle, sigma(angle)
- GSASIImath.getAngleDerv(Oxyz, Axyz, Bxyz, Amat, Tunit, symNo, SGData)[source]
computes the numerical derivative of the angle between two vectors (generated from pairs of atoms) used here in
getAngSig()
.
- GSASIImath.getAtomPtrs(data, draw=False)[source]
get atom data pointers cx,ct,cs,cia in Atoms or Draw Atoms lists NB:may not match column numbers in displayed table
param: dict: data phase data structure draw: boolean True if Draw Atoms list pointers are required return: cx,ct,cs,cia pointers to atom xyz, type, site sym, uiso/aniso flag
- GSASIImath.getAtomXYZ(atoms, cx)[source]
Create an array of fractional coordinates from the atoms list
- Parameters:
atoms (list) – atoms object as found in tree
cx (int) – offset to where coordinates are found
- Returns:
np.array with shape (n,3)
- GSASIImath.getCWgam(ins, pos)[source]
get CW peak profile gamma
- Parameters:
ins (dict) – instrument parameters with at least ‘X’, ‘Y’ & ‘Z’ as values only
pos (float) – 2-theta of peak
- Returns:
float getCWgam: peak gamma
- GSASIImath.getCWgamDeriv(pos)[source]
get derivatives of CW peak profile gamma wrt X, Y & Z
- Parameters:
pos (float) – 2-theta of peak
- Returns:
list getCWgamDeriv: d(gam)/dX & d(gam)/dY
- GSASIImath.getCWsig(ins, pos)[source]
get CW peak profile sigma^2
- Parameters:
ins (dict) – instrument parameters with at least ‘U’, ‘V’, & ‘W’ as values only
pos (float) – 2-theta of peak
- Returns:
float getCWsig: peak sigma^2
- GSASIImath.getCWsigDeriv(pos)[source]
get derivatives of CW peak profile sigma^2 wrt U,V, & W
- Parameters:
pos (float) – 2-theta of peak
- Returns:
list getCWsigDeriv: d(sig^2)/dU, d(sig)/dV & d(sig)/dW
- GSASIImath.getDensity(generalData)[source]
calculate crystal structure density
- Parameters:
generalData (dict) – The General dictionary in Phase
- Returns:
float density: crystal density in gm/cm^3
- GSASIImath.getDistDerv(Oxyz, Txyz, Amat, Tunit, Top, SGData)[source]
computes the numerical derivative of the distance between two atoms used here in
CalcDistDeriv()
and inGSASIIstrMain.RetDistAngle()
- Parameters:
name (type) – description
- Returns:
type name: description
- GSASIImath.getEDgam(ins, pos)[source]
get ED peak profile gam
- Parameters:
ins (dict) – instrument parameters with at least X, Y & Z as values only
pos (float) – energy of peak as keV
- Returns:
float getEDsig: peak gam im keV
- GSASIImath.getEDgamDeriv(ins, pos)[source]
get derivatives of ED peak profile gam wrt X, Y & Z
- Parameters:
pos (float) – energy of peak in keV
- Returns:
list getEDsigDeriv: d(gam)/dX, d(gam)/dY & d(gam)/dZ,
- GSASIImath.getEDsig(ins, pos)[source]
get ED peak profile sig
- Parameters:
ins (dict) – instrument parameters with at least ‘A’, ‘B’ & ‘C’ as values only
pos (float) – energy of peak as keV
- Returns:
float getEDsig: peak sigma^2 im keV**2
- GSASIImath.getEDsigDeriv(ins, pos)[source]
get derivatives of ED peak profile sig wrt A, B & C
- Parameters:
pos (float) – energy of peak in keV
- Returns:
list getEDsigDeriv: d(sig)/dA, d(sig)/dB & d(sig)/dC,
- GSASIImath.getMass(generalData)[source]
Computes mass of unit cell contents
- Parameters:
generalData (dict) – The General dictionary in Phase
- Returns:
float mass: Crystal unit cell mass in AMU.
- GSASIImath.getMeanWave(Parms)[source]
returns mean wavelength from Instrument parameters dictionary
- Parameters:
Parms (dict) – Instrument parameters; must contain: Lam: single wavelength or Lam1,Lam2: Ka1,Ka2 radiation wavelength I(L2)/I(L1): Ka2/Ka1 ratio
- Returns:
float wave: mean wavelength
- GSASIImath.getPinkAlpha(ins, tth)[source]
get pink neutron peak alpha profile
- Parameters:
ins (dict) – instrument parameters with at least ‘alpha’ as values only
tth (float) – 2-theta of peak
- Returns:
float getPinkNalpha: peak alpha
- GSASIImath.getPinkAlphaDeriv(tth)[source]
get alpha derivatives of pink neutron peak profile
- Parameters:
tth (float) – 2-theta of peak
- Returns:
float getPinkNalphaDeriv: d(alp)/d(alpha-0), d(alp)/d(alpha-1)
- GSASIImath.getPinkBeta(ins, tth)[source]
get pink neutron peak profile beta
- Parameters:
ins (dict) – instrument parameters with at least ‘beta-0’ & ‘beta-1’ as values only
tth (float) – 2-theta of peak
- Returns:
float getPinkbeta: peak beta
- GSASIImath.getPinkBetaDeriv(tth)[source]
get beta derivatives of pink neutron peak profile
- Parameters:
tth (float) – 2-theta of peak
- Returns:
list getPinkNbetaDeriv: d(beta)/d(beta-0) & d(beta)/d(beta-1)
- GSASIImath.getRBTransMat(X, Y)[source]
Get transformation for Cartesian axes given 2 vectors X will be parallel to new X-axis; X cross Y will be new Z-axis & (X cross Y) cross Y will be new Y-axis Useful for rigid body axes definintion
- Parameters:
X (array) – normalized vector
Y (array) – normalized vector
- Returns:
array M: transformation matrix
use as XYZ’ = np.inner(M,XYZ) where XYZ are Cartesian
- GSASIImath.getRamaDeriv(XYZ, Amat, Coeff)[source]
Computes numerical derivatives of torsion angle pair pseudo potential with respect of crystallographic atom coordinates of the 5 atom sequence
- Parameters:
XYZ (nparray) – crystallographic coordinates of 5 atoms
Amat (nparray) – crystal to cartesian transformation matrix
Coeff (list) – pseudo potential coefficients
- Returns:
list (deriv) derivatives of pseudopotential with respect to 5 atom crystallographic xyz coordinates.
- GSASIImath.getRestAngle(XYZ, Amat)[source]
Compute interatomic angle(s) for use in restraints
- Parameters:
name (type) – description
- Returns:
type name: description
- GSASIImath.getRestChiral(XYZ, Amat)[source]
compute a chiral restraint
- Parameters:
name (type) – description
- Returns:
type name: description
- GSASIImath.getRestDeriv(Func, XYZ, Amat, ops, SGData)[source]
Compute numerical derivatives of restraints for use in minimization
- Parameters:
name (type) – description
- Returns:
type name: description
- GSASIImath.getRestDist(XYZ, Amat)[source]
Compute interatomic distance(s) for use in restraints
- Parameters:
name (type) – description
- Returns:
type name: description
- GSASIImath.getRestPlane(XYZ, Amat)[source]
Compute deviations from a best plane through atoms for use in restraints
- Parameters:
name (type) – description
- Returns:
type name: description
- GSASIImath.getRestPolefig(ODFln, SamSym, Grid)[source]
default doc string
- Parameters:
name (type) – description
- Returns:
type name: description
- GSASIImath.getRestPolefigDerv(HKL, Grid, SHCoeff)[source]
default doc string
- Parameters:
name (type) – description
- Returns:
type name: description
- GSASIImath.getRestRama(XYZ, Amat)[source]
Computes a pair of torsion angles in a 5 atom string
- Parameters:
XYZ (nparray) – crystallographic coordinates of 5 atoms
Amat (nparray) – crystal to cartesian transformation matrix
- Returns:
list (phi,psi) two torsion angles in degrees
- GSASIImath.getRestTorsion(XYZ, Amat)[source]
compute a torsion restraint
- Parameters:
name (type) – description
- Returns:
type name: description
- GSASIImath.getRho(xyz, mapData)[source]
get scattering density at a point by 8-point interpolation param xyz: coordinate to be probed param: mapData: dict of map data
- Returns:
density at xyz
- GSASIImath.getRhos(XYZ, rho)[source]
get scattering density at an array of point by 8-point interpolation this is faster than gerRho which is only used for single points. However, getRhos is replaced by scipy.ndimage.interpolation.map_coordinates which does a better job & is just as fast. Thus, getRhos is unused in GSAS-II at this time. param xyz: array coordinates to be probed Nx3 param: rho: array copy of map (NB: don’t use original!)
- Returns:
density at xyz
- GSASIImath.getSyXYZ(XYZ, ops, SGData)[source]
default doc
- Parameters:
name (type) – description
- Returns:
type name: description
- GSASIImath.getTOFalpha(ins, dsp)[source]
get TOF peak profile alpha
- Parameters:
ins (dict) – instrument parameters with at least ‘alpha’ as values only
dsp (float) – d-spacing of peak
- Returns:
float getTOFalpha: peak alpha
- GSASIImath.getTOFalphaDeriv(dsp)[source]
get alpha derivatives of TOF peak profile
- Parameters:
dsp (float) – d-spacing of peak
- Returns:
float getTOFalphaDeriv: d(alp)/d(alpha)
- GSASIImath.getTOFbeta(ins, dsp)[source]
get TOF peak profile beta
- Parameters:
ins (dict) – instrument parameters with at least ‘beat-0’, ‘beta-1’ & ‘beta-q’ as values only
dsp (float) – d-spacing of peak
- Returns:
float getTOFbeta: peak beat
- GSASIImath.getTOFbetaDeriv(dsp)[source]
get derivatives of TOF peak profile beta wrt beta-0, beta-1, & beat-q
- Parameters:
dsp (float) – d-spacing of peak
- Returns:
list getTOFbetaDeriv: d(beta)/d(beat-0), d(beta)/d(beta-1) & d(beta)/d(beta-q)
- GSASIImath.getTOFgamma(ins, dsp)[source]
get TOF peak profile gamma
- Parameters:
ins (dict) – instrument parameters with at least ‘X’, ‘Y’ & ‘Z’ as values only
dsp (float) – d-spacing of peak
- Returns:
float getTOFgamma: peak gamma
- GSASIImath.getTOFgammaDeriv(dsp)[source]
get derivatives of TOF peak profile gamma wrt X, Y & Z
- Parameters:
dsp (float) – d-spacing of peak
- Returns:
list getTOFgammaDeriv: d(gam)/dX & d(gam)/dY
- GSASIImath.getTOFsig(ins, dsp)[source]
get TOF peak profile sigma^2
- Parameters:
ins (dict) – instrument parameters with at least ‘sig-0’, ‘sig-1’ & ‘sig-q’ as values only
dsp (float) – d-spacing of peak
- Returns:
float getTOFsig: peak sigma^2
- GSASIImath.getTOFsigDeriv(dsp)[source]
get derivatives of TOF peak profile sigma^2 wrt sig-0, sig-1, & sig-q
- Parameters:
dsp (float) – d-spacing of peak
- Returns:
list getTOFsigDeriv: d(sig0/d(sig-0), d(sig)/d(sig-1) & d(sig)/d(sig-q)
- GSASIImath.getTorsionDeriv(XYZ, Amat, Coeff)[source]
default doc string
- Parameters:
name (type) – description
- Returns:
type name: description
- GSASIImath.getVCov(varyNames, varyList, covMatrix)[source]
obtain variance-covariance terms for a set of variables. NB: the varyList and covMatrix were saved by the last least squares refinement so they must match.
- Parameters:
varyNames (list) – variable names to find v-cov matric for
varyList (list) – full list of all variables in v-cov matrix
covMatrix (nparray) – full variance-covariance matrix from the last least squares refinement
- Returns:
nparray vcov: variance-covariance matrix for the variables given in varyNames
- GSASIImath.getWave(Parms)[source]
returns wavelength from Instrument parameters dictionary
- Parameters:
Parms (dict) – Instrument parameters; must contain: Lam: single wavelength or Lam1: Ka1 radiation wavelength
- Returns:
float wave: wavelength
- GSASIImath.make2Quat(A, B)[source]
Make quaternion from rotation of A vector to B vector
- Parameters:
A,B (np.array) – Cartesian 3-vectors
- Returns:
quaternion & rotation angle in radians q=r+ai+bj+ck
- GSASIImath.makeQuat(A, B, C)[source]
Make quaternion from rotation of A vector to B vector about C axis
- Parameters:
A,B,C (np.array) – Cartesian 3-vectors
- Returns:
quaternion & rotation angle in radians q=r+ai+bj+ck
- GSASIImath.makeWaves(waveTypes, FSSdata, XSSdata, USSdata, MSSdata, Mast)[source]
waveTypes: array nAtoms: ‘Fourier’,’ZigZag’ or ‘Block’ FSSdata: array 2 x atoms x waves (sin,cos terms) XSSdata: array 2x3 x atoms X waves (sin,cos terms) USSdata: array 2x6 x atoms X waves (sin,cos terms) MSSdata: array 2x3 x atoms X waves (sin,cos terms)
Mast: array orthogonalization matrix for Uij
- GSASIImath.makeWavesDerv(ngl, waveTypes, FSSdata, XSSdata, USSdata, Mast)[source]
Only for Fourier waves for fraction, position & adp (probably not used for magnetism) FSSdata: array 2 x atoms x waves (sin,cos terms) XSSdata: array 2x3 x atoms X waves (sin,cos terms) USSdata: array 2x6 x atoms X waves (sin,cos terms) Mast: array orthogonalization matrix for Uij
- GSASIImath.mcsaSearch(data, RBdata, reflType, reflData, covData, pgbar, start=True)[source]
default doc string
- Parameters:
name (type) – description
- Returns:
type name: description
- GSASIImath.phaseContents(phase)[source]
Compute the unit cell and asymmetric unit contents for a phase.
This has been tested only on type=’nuclear’ phases and might need adaptation for phases of other types, if the phase type does not have an occupancy defined.
- Parameters:
phase (dict) – the dict for a phase, as found in the data tree
- Returns:
acomp,ccomp where acomp is the asymmetric unit contents and ccomp is the contents of the unit cell
- GSASIImath.pinv(a, rcond=1e-15)[source]
Compute the (Moore-Penrose) pseudo-inverse of a matrix. Modified from numpy.linalg.pinv; assumes a is Hessian & returns no. zeros found Calculate the generalized inverse of a matrix using its singular-value decomposition (SVD) and including all large singular values.
- Parameters:
a (array) – (M, M) array_like - here assumed to be LS Hessian Matrix to be pseudo-inverted.
rcond (float) – Cutoff for small singular values. Singular values smaller (in modulus) than rcond * largest_singular_value (again, in modulus) are set to zero.
- Returns:
B : (M, M) ndarray The pseudo-inverse of a
- Raises: LinAlgError
If the SVD computation does not converge.
- Notes:
The pseudo-inverse of a matrix A, denoted \(A^+\), is defined as: “the matrix that ‘solves’ [the least-squares problem] \(Ax = b\),” i.e., if \(\bar{x}\) is said solution, then \(A^+\) is that matrix such that \(\bar{x} = A^+b\).
It can be shown that if \(Q_1 \Sigma Q_2^T = A\) is the singular value decomposition of A, then \(A^+ = Q_2 \Sigma^+ Q_1^T\), where \(Q_{1,2}\) are orthogonal matrices, \(\Sigma\) is a diagonal matrix consisting of A’s so-called singular values, (followed, typically, by zeros), and then \(\Sigma^+\) is simply the diagonal matrix consisting of the reciprocals of A’s singular values (again, followed by zeros). [1]
References: .. [1] G. Strang, Linear Algebra and Its Applications, 2nd Ed., Orlando, FL, Academic Press, Inc., 1980, pp. 139-142.
- GSASIImath.printRho(SGLaue, rho, rhoMax)[source]
default doc string
- Parameters:
name (type) – description
- Returns:
type name: description
- GSASIImath.randomAVdeg(r0, r1, r2, r3)[source]
create random angle (deg),vector from 4 random number in range (-1,1)
- GSASIImath.randomQ(r0, r1, r2, r3)[source]
create random quaternion from 4 random numbers in range (-1,1)
- GSASIImath.searchBondRestr(origAtoms, targAtoms, bond, Factor, GType, SGData, Amat, defESD=0.01, dlg=None)[source]
Search for bond distance restraints.
- GSASIImath.setHcorr(info, Amat, xtol, problem=False)[source]
Find & report high correlation terms in covariance matrix
- GSASIImath.setPeakparms(Parms, Parms2, pos, mag, ifQ=False, useFit=False)[source]
set starting peak parameters for single peak fits from plot selection or auto selection
- Parameters:
Parms (dict) – instrument parameters dictionary
Parms2 (dict) – table lookup for TOF profile coefficients
pos (float) – peak position in 2-theta, TOF or Q (ifQ=True)
mag (float) – peak top magnitude from pick
ifQ (bool) – True if pos in Q
useFit (bool) – True if use fitted CW Parms values (not defaults)
- Returns:
list XY: peak list entry: for CW: [pos,0,mag,1,sig,0,gam,0] for TOF: [pos,0,mag,1,alp,0,bet,0,sig,0,gam,0] for Pink: [pos,0,mag,1,alp,0,bet,0,sig,0,gam,0] NB: mag refinement set by default, all others off