57Fe Mossbauer spectrum of DO3 alpha-Fe-Si phase

Fit function libraries that can be linked to and utilized in the MossWinn program.
Z. Klencsar
Posts: 9
Joined: Fri May 20, 2016 11:17 pm

57Fe Mossbauer spectrum of DO3 alpha-Fe-Si phase

Postby Z. Klencsar » Tue May 24, 2016 5:25 pm

library sub_finemetdo3a;

{

Written by Zoltán Klencsár, z.klencsar@mosswinn.hu

The routine calculates contribution of the nanocrystalline alpha-Fe-Si phase to the 57Fe Mossbauer spectrum
of annealed FINEMET as a function of Si concentration of the DO3 alpha-Fe-Si phase.

The code can be compiled with 32 bit Delphi compilers to a DLL library that can be used together
with the MossWinn program (version 4.0Pre or later) to fit Mossbauer spectra.

The source code of this library is also available for download from
http://www.mosswinn.hu/dlls/sub_finemetdo3a.dpr

A compiled version of this code is available for download from
http://www.mosswinn.hu/dlls/sub_finemetdo3a.dll

Copy the compiled DLL into the folder ..\MossWinn 4.0\DLLs\ in order to be able to use it in MossWinn.

If you use this code to fit spectra and publish the results, please, in your work reference the corresponding paper of

Oshtrakh et al.: Materials Chemistry and Physics (2016).

http://dx.doi.org/10.1016/j.matchemphys.2016.05.032

If you improve upon this library code, please, post the improvement here as answer.

}

uses Math;

procedure DLL1_SET_INTERACTION_NAME(var s:ShortString;source:ShortString);export;
begin
if (source='57FE') then s:='Fe-Si DO3 MODEL A' else s:='Undefined';
end;

procedure DLL1_SET_DESCRIPTION(var s:ShortString;source:ShortString);export;
begin
DLL1_SET_INTERACTION_NAME(s,source);{For simplicity description is set equal to the interaction name.}
end;

{Setting number of parameters}
procedure DLL1_SET_SUBSPECTRUM_INTERNAL_VALUES(var AN,PN,WN:byte;
var Name:ShortString;Source:ShortString;
p:pointer);export;
begin
AN:=2;{Amplitude type}
PN:=17;{Position type}
WN:=6;{Width type}
if (source='57FE') then Name:='Fe-Si DO3 MODEL A (1)' else Name:='Undefined (1)';
end;

{Initialization of parameters}
procedure DLL1_SET_PARAMETERS_INTERNAL_VALUES(n:longint;
Cmin,Cmax,Vmin,Vmax:double;
var Name:ShortString;
var Minimum_Value,
Maximum_Value,
Initial_Value:double;
Source:ShortString;
p:pointer);export;
begin

Case n of
1:begin
Name:='TOTAL AREA (DO3)';{Total area of the component associated with alpha-Fe-Si.}
Minimum_Value:=0;
Maximum_Value:=(Vmax-Vmin)*(Cmax-Cmin)/4;
Initial_Value:=Maximum_Value/5;
end;

2:begin
Name:='A25 / A34';{Relative area of the 2-5 peaks of the sextets with respect to that of the 3-4 peaks.}
Minimum_Value:=0;
Maximum_Value:=4;
Initial_Value:=4;
end;

3:begin
Name:='ISOMER SHIFT (A4)'; {A site iron with 4 n.n. D site Fe atoms.}
Minimum_Value:=0.22;
Maximum_Value:=0.28;
Initial_Value:=0.25;
end;

4:begin
Name:='ISOMER SHIFT (A5)'; {A site iron with 5 n.n. D site Fe atoms.}
Minimum_Value:=0.17;
Maximum_Value:=0.23;
Initial_Value:=0.18;
end;

5:begin
Name:='ISOMER SHIFT (A6)'; {A site iron with 6 n.n. D site Fe atoms.}
Minimum_Value:=-0.05;
Maximum_Value:=0.15;
Initial_Value:=0.1;
end;

6:begin
Name:='ISOMER SHIFT (A7)'; {A site iron with 7 n.n. D site Fe atoms.}
Minimum_Value:=-0.05;
Maximum_Value:=0.12;
Initial_Value:=0.05;
end;

7:begin
Name:='ISOMER SHIFT (D0)';
Minimum_Value:=-0.05;
Maximum_Value:=0.15;
Initial_Value:=0.07;
end;

8:begin
Name:='ISOMER SHIFT (Dprime+A8)';{D'+A8}
Minimum_Value:=-0.05;
Maximum_Value:=0.15;
Initial_Value:=0.03;
end;

9:begin
Name:='HYP. M. FIELD (A4)';
Minimum_Value:=18;
Maximum_Value:=22;
Initial_Value:=19.7;
end;

10:begin
Name:='HYP. M. FIELD (A5)';
Minimum_Value:=22;
Maximum_Value:=27;
Initial_Value:=24.6;
end;

11:begin
Name:='HYP. M. FIELD (A6)';
Minimum_Value:=27;
Maximum_Value:=30;
Initial_Value:=28.9;
end;

12:begin
Name:='HYP. M. FIELD (A7)';
Minimum_Value:=30;
Maximum_Value:=31.5;
Initial_Value:=30.2;
end;

13:begin
Name:='HYP. M. FIELD (D0)';
Minimum_Value:=31.0;
Maximum_Value:=32.0;
Initial_Value:=31.1;
end;

14:begin
Name:='HYP. M. FIELD (Dprime+A8)';
Minimum_Value:=32.0;
Maximum_Value:=33.0;
Initial_Value:=32.4;
end;

15:begin
Name:='2*EPSILON (A4)'; {2*quadrupole shift}
Minimum_Value:=-0.1;
Maximum_Value:=0.1;
Initial_Value:=0.0;
end;

16:begin
Name:='2*EPSILON (A5)';
Minimum_Value:=-0.1;
Maximum_Value:=0.1;
Initial_Value:=0.0;
end;

17:begin
Name:='2*EPSILON (A6)';
Minimum_Value:=-0.1;
Maximum_Value:=0.1;
Initial_Value:=0.0;
end;

18:begin
Name:='2*EPSILON (A7)';
Minimum_Value:=-0.1;
Maximum_Value:=0.1;
Initial_Value:=0.0;
end;

19:begin
Name:='Si concentration (DO3)';
Minimum_Value:=15;
Maximum_Value:=25;
Initial_Value:=20.5;
end;

20:begin
Name:='LINE WIDTH (A4)';
Minimum_Value:=0.23;
Maximum_Value:=0.40;
Initial_Value:=0.31;
end;

21:begin
Name:='LINE WIDTH (A5)';
Minimum_Value:=0.23;
Maximum_Value:=0.40;
Initial_Value:=0.31;
end;

22:begin
Name:='LINE WIDTH (A6)';
Minimum_Value:=0.23;
Maximum_Value:=0.40;
Initial_Value:=0.27;
end;

23:begin
Name:='LINE WIDTH (A7)';
Minimum_Value:=0.23;
Maximum_Value:=0.40;
Initial_Value:=0.27;
end;

24:begin
Name:='LINE WIDTH (D0)';
Minimum_Value:=0.23;
Maximum_Value:=0.40;
Initial_Value:=0.35;
end;

25:begin
Name:='LINE WIDTH (Dprime+A8)';
Minimum_Value:=0.23;
Maximum_Value:=0.40;
Initial_Value:=0.31;
end;

end;

end;

procedure DLL1_SET_SUBSPECTRUM(n:longint;
pPARAMETERS,pDATA,pVELOCITY:pointer;
AreaMode:boolean;Source:ShortString;p:pointer);export;
CONST
twopi = pi+pi;

{Physical constants}
speed_of_light = 2.99792458E8;
elementary_charge = 1.60217653E-19;

gFe32 = -0.103267; {I=3/2 giromagnetic factor of 57Fe}
gFe12 = +0.18088; {I=1/2 giromagnetic factor of 57Fe}

nuclear_magneton = 5.05078343E-27; {Joule/Tesla}
hreduced = 1.05457168E-34; {Js/rad} {Planck constant divided by 2*pi}

EgammaFe = 14.4129; {keV}
Fe_Egamma = EgammaFe*1000.0*elementary_charge; {Joule}
mmpstoJoule57Fe = 0.001*Fe_Egamma/speed_of_light;

TYPE
TARRAY = array[1..8192] of double;

{Mossbauer parameters characteristic of a Fe microenvironment.}
TSiteRecord = Record
Weight:double;
ISH:double;
QS:double;
HMF:double;
LW:double;
End;

VAR
HMe,HMg:double;
PARAMETERS,DATA,VELOCITY:^TARRAY;
Area,RELA:double;
A4,A5,A6,A7,D0,Dprime:TSiteRecord;
CSi:double;
PB1,PB2,PB3,PB4,PB5,PB6:double;
WeightSUM:double;

{Adds a Lorentzian singlet with area A, center V0 and FWHM LW.}
Procedure AddSinglet(const A,V0,LW:double);
var
i:longint;
xa,x,y,v,LW2p4:double;
begin
xa:=A*LW/twopi;
LW2p4:=LW*LW/4;
For i:=1 to n do
begin
v:=VELOCITY^[i];
x:=(v-v0);
y:=x*x+LW2p4;
DATA^[i]:=DATA^[i]+xa/y;
end;
end;

{Adds a sextet with full area A, isomer shift ISH, hyperfine magnetic field HMF,
2*quadrupole shift of TWO_EPSILON, FWHM of LW.}
Procedure AddSextet(const A,ISH,HMF,TWO_EPSILON,LW:double);
var
EPSILON:double;
LP1,LP2,LP3,LP4,LP5,LP6:double;
LA16,LA25,LA34:double;
begin
EPSILON:=TWO_EPSILON/2.0;

{Calculation of the position of individual peaks of the sextet.}
LP1:=ISH+EPSILON+PB1*HMF;
LP2:=ISH-EPSILON+PB2*HMF;
LP3:=ISH-EPSILON+PB3*HMF;
LP4:=ISH-EPSILON+PB4*HMF;
LP5:=ISH-EPSILON+PB5*HMF;
LP6:=ISH+EPSILON+PB6*HMF;

{Calculation of individual peak areas from the total area of the sextet and A25/A34.}
LA34:=A/(8+RELA+RELA);
LA16:=3*LA34;
LA25:=RELA*LA34;

{Adding the 6 peaks of the sextet.}
AddSinglet(LA16,LP1,LW);AddSinglet(LA16,LP6,LW);
AddSinglet(LA25,LP2,LW);AddSinglet(LA25,LP5,LW);
AddSinglet(LA34,LP3,LW);AddSinglet(LA34,LP4,LW);
end;

Function WeightOfSite(const SiteNumber:longint):double;
var
SiPercent:double;
psiteD:double;
binomialfactor,xx:double;
x,pDprime_per_pd0,pd0:double;
begin
SiPercent:=CSi;
if SiPercent>24.99999 then SiPercent:=24.99999; {To avoid Run Time Error for ln(0)}
{The code is not valid for SiPercent>=25.}

case SiteNumber of
4:binomialfactor:=1;{A4}
5:binomialfactor:=4;{A5}
6:binomialfactor:=6;{A6}
7:binomialfactor:=4;{A7}
8:binomialfactor:=1;{A8}
end;

{Fe2 [ Fe1+x Si1-x ]} {SiPercent = 100*(1-x)/4}
xx:=SiPercent*(2/100); { = (1-x)/2 }
pSiteD:=1.0-xx; { = (1+x)/2 } {Probability of occurrence of iron at D sites.}

x:=1.0-xx-xx;
pDprime_per_pd0:=(1+x)/((1-x)*(1-x)*(1-x)*(1-x)*(1-x)*(1+5*x))-1;{See Eqs. (7) and (9).}
pd0:=pSiteD/(pDprime_per_pd0+1);

if SiteNumber=9 then RESULT:=pDprime_per_pd0 {This gives relative area of Fe at D sites of Si + D sites of Fe where at neighboring Si sites there are 2 or more iron atoms. This will have higher hyperfine magnetic field than that of iron at D sites with 1 or 0 D site Fe neighbor.}
else RESULT:=binomialfactor*EXP((8-SiteNumber)*ln(1.0-x))*EXP((SiteNumber-4)*ln(x))/pd0; {See Eq.(1)}
end;

begin
if NOT(source='57FE') then exit;

HMe:=-gFe32*nuclear_magneton/mmpstoJoule57Fe;
HMg:=-gFe12*nuclear_magneton/mmpstoJoule57Fe;

{Sextet line positions for 1 T hyperfine magnetic field, in mm/s}
PB1:=HMe*(-3/2)-HMg*(-1/2);
PB2:=HMe*(-1/2)-HMg*(-1/2);
PB3:=HMe*(1/2)-HMg*(-1/2);
PB4:=HMe*(-1/2)-HMg*(1/2);
PB5:=HMe*(1/2)-HMg*(1/2);
PB6:=HMe*(3/2)-HMg*(1/2);

PARAMETERS:=pPARAMETERS;DATA:=pDATA;VELOCITY:=pVELOCITY;

Area:=MAX(PARAMETERS^[1],0);{Total Area}
RELA:=MIN(MAX(PARAMETERS^[2],0),4);{Relative area of 2-5 peaks with respect to the 3-4 peaks for all sextets.}

WITH A4 DO begin ISH:=PARAMETERS^[3];HMF:=PARAMETERS^[9];QS:=PARAMETERS^[15];LW:=MAX(PARAMETERS^[20],0.01); end;
WITH A5 DO begin ISH:=PARAMETERS^[4];HMF:=PARAMETERS^[10];QS:=PARAMETERS^[16];LW:=MAX(PARAMETERS^[21],0.01); end;
WITH A6 DO begin ISH:=PARAMETERS^[5];HMF:=PARAMETERS^[11];QS:=PARAMETERS^[17];LW:=MAX(PARAMETERS^[22],0.01); end;
WITH A7 DO begin ISH:=PARAMETERS^[6];HMF:=PARAMETERS^[12];QS:=PARAMETERS^[18];LW:=MAX(PARAMETERS^[23],0.01); end;
WITH D0 DO begin ISH:=PARAMETERS^[7];HMF:=PARAMETERS^[13];QS:=0.0;LW:=MAX(PARAMETERS^[24],0.01); end; {Fe at regular D Fe sites with 0 or 1 n.n.n Fe at D Si sites.}
WITH Dprime DO begin ISH:=PARAMETERS^[8];HMF:=PARAMETERS^[14];QS:=0.0;LW:=MAX(PARAMETERS^[25],0.01); end; {Fe at Si D sites + A8 + Fe at regular D Fe sites with 2 to 6 n.n.n Fe at D Si sites.}

CSi:=PARAMETERS^[19];{Si concentration}

A4.Weight:=WeightOfSite(4);
A5.Weight:=WeightOfSite(5);
A6.Weight:=WeightOfSite(6);
A7.Weight:=WeightOfSite(7);
Dprime.Weight:=WeightOfSite(8){A8}+WeightOfSite(9);
D0.Weight:=1.0;{The previous ones are relative to this value.}

WeightSUM:=A4.Weight+A5.Weight+A6.Weight+A7.Weight+Dprime.Weight+D0.Weight;

{Now we make the sum of weights to be equal to 1.0.}
A4.Weight:=A4.Weight/WeightSUM;
A5.Weight:=A5.Weight/WeightSUM;
A6.Weight:=A6.Weight/WeightSUM;
A7.Weight:=A7.Weight/WeightSUM;
Dprime.Weight:=Dprime.Weight/WeightSUM;
D0.Weight:=D0.Weight/WeightSUM;

With A4 do AddSextet(Weight*Area,ISH,HMF,QS,LW);
With A5 do AddSextet(Weight*Area,ISH,HMF,QS,LW);
With A6 do AddSextet(Weight*Area,ISH,HMF,QS,LW);
With A7 do AddSextet(Weight*Area,ISH,HMF,QS,LW);
With D0 do AddSextet(Weight*Area,ISH,HMF,QS,LW);
With Dprime do AddSextet(Weight*Area,ISH,HMF,QS,LW);
end;


{Set spectrum area.}
procedure DLL1_SET_AREA_OF_SUBSPECTRUM(var AREA:double;
var CanDO:boolean;
pPARAMETERS:pointer;
AreaMode:boolean;Source:ShortString;
p:pointer);export;
TYPE
TARRAY = array[1..8192] of double;
var
PARAMETERS:^TARRAY;
begin
CanDO:=True;

if (source='57FE') then
begin
PARAMETERS:=pPARAMETERS;
Area:=MAX(PARAMETERS^[1],0);
end
else Area:=0.0;
end;

exports

DLL1_SET_INTERACTION_NAME,
DLL1_SET_DESCRIPTION,
DLL1_SET_SUBSPECTRUM_INTERNAL_VALUES,
DLL1_SET_PARAMETERS_INTERNAL_VALUES,
DLL1_SET_SUBSPECTRUM,
DLL1_SET_AREA_OF_SUBSPECTRUM;

begin

end.