ross-colman
Hi,
I'm possibly trying to do more than is reasonable, but anyway:
I'm attempting the exhaustive symmetry analysis approach to fitting some of my own data, but I'm including site occupancy variation too (checking for cation ordering of a site that is mixed occupancy (disordered) in the parent structure.
As I'm allowing the occupancy to vary, and as the occupancy is a dependent parameter (with the symmetry modes being the primary parameters that determine it), I can't put a min-max into each site line. In some of the structures, it readily goes to unphysical values.
I'd like to put in a penalty function that pushes the total occupancies within the unit cell towards the ideal formula.
As the exhaustive symmetry analysis runs through multiple structure files, with each structure having a different number of sites and different unit cell parameters, I'm struggling to figure the best way to do it.
Ideally it will be a macro within the starting_parent.inp file.
I think I will have to set up a different penalty function for each atom type as there will be different numbers of sites for different atoms within the various .str files.
My idea was to add up the total number of atoms in the unit cell, then divide by the volume to get the 'species density': species density = (sum over all sites of the atom (num_posns * occ)) / cell volume
Then penalize variation away from the ideal species density.
My particular difficulty is in doing this with varying numbers of sites of the specific atoms - so need some sort of loop.
I took inspiration from the 'cif file writing macro' and made a first attempt at a solution but unsurprisingly it doesn't work (i suppose because these commands were designed for writing to files, not for use in an equation) and I'm not sure what to try now:
"
macro cP(species, ideal_species_density)
{
atom_out
{
local species_density = 0;: 0
existing_prm species_density += IF(Get_From_String(Get(current_atom), atom) == species , (Get_From_String(Get(current_atom), num_posns) * Get_From_String(Get(current_atom), occ)) / Get(V) ; , 0);
local species_density_difference = ideal_species_density - species_density;: 0
}
penalty = Ln(Abs(species_density_difference)+1);
}
"
In my defense, this is possibly the first topas macro I've ever attempted to write.
Has anyone got any advice on how to approach this problem?
Thinking about it now, it might be easier to get the python script to construct the penalty whilst translating each .str file into a .inp - but in case there is a quick way to get this macro idea to work then I'd love to hear it.
Kind regards,
Ross
rowlesmr
Why not just calculate the density from the unit cell mass and volume?
This is an already existing macro in topas.inc:
macro & Phase_Density_Eqn_g_on_cm3
{
1.6605402 Get(cell_mass) / Get(cell_volume)
}
Matthew
* 1.6605402 = 1/N_A to convert Da/Å^3 to g/cm^3.
ross-colman
Dear Matthew,
Thank you for the suggestion.
The difficulty is that simply controlling the density won't necessarily constrain the composition enough.
There is the potential for both A3+/A4+ order/disorder, as well as oxygen/vacancy order/disorder within the subgroups that I'm trialing. They could cross-compensate, so I think I will need to include atom-specific constraints.
My issue is still on how to sum over all sites of one atom type, when different subgroups will have different numbers of sites for that atom.
I'd welcome any other suggestions or thoughts you had on how to do it though.
Many thanks,
Ross
johnsoevans
Ross,
Do any of the overall composition commands like element_weight_percent do what you want? If you're using jEdit look in the Examples/Quantitative folder for examples. These let you add up overall weight percentages across multiple phases so normally go at the xdd level. You can then switch to composition with mass/RMM type equations.
I've also pasted a few lines from an INP file below. Experience suggests you'll have to test refinement stability carefully if using these.
John
element_weight_percent Sr wpcSr 13.82175404`_0.35551868
element_weight_percent Ca wpcCa 13.84611784`_0.38670822
element_weight_percent Nb wpcNb 48.02592037`_0.28310688
element_weight_percent O wpcO 24.30620775`_0.13683736
prm mol_Sr =wpcSr/87.620;:0.157746565`_0.00405750606
prm mol_Ca =wpcCa/137.33;:0.100823694`_0.00281590491
prm mol_Nb =wpcNb/92.906;:0.516930235`_0.00304724005
prm mol_O =wpcO /15.999;:1.51923294`_0.00855286977
prm comp_Sr = (mol_Sr/mol_Nb)*2;:0.610320521`_0.0177510082
prm comp_Ca = (mol_Ca/mol_Nb)*2;:0.39008627`_0.0124681274
prm comp_Nb = (mol_Nb/mol_Nb)*2;:2`_0
prm comp_O = (mol_O /mol_Nb)*2;:5.87790319`_0.00560845677
'Add penalties here to control the overall SrxCa1-x in all phases
penalties_weighting_K1 1
penalty = (70 - 100 * comp_Sr)^2;:0.00102734026`
penalty = (30 - 100 * comp_Ca)^2;:7.44255972e-05`
prm test = ((70/30) - (comp_Sr/comp_Ca))^2;:2.26343906e-07`_6.1787783e-05
ross-colman
Dear John,
Thank you, I've just done a very quick test and I think this is exactly what I need!
It should allow me to constrain the overall composition irrespective of how many sites there are.
You have my gratitude again!
Ross