1. Batch ThermoCalc Executables

In the present study we are interested in seeing how nitrogen, and titanium affect the resulting microstructure, as well as their interactions with the other components in a 2032Nb alloy. From the ASTM specification for 20Cr-32Nb-Nb the compositional matrix in Table 1 can be used see the ranges for each element that will be explored. In ThermoCalc the stepping command will be used over the temperature range T_m - T_rt for finite compositions of the alloy. Therefore, certain levels in the compositional ranges for each element will need to be decided where then equilibrium will be calculated for all possible combinations of these finite chemistries. If a maximum, median, and minimum values are chosen for each range than ThermoCalc will need to run stepping procedures for 3⁷ = 2187 chemical compositions. Writing these ThermoCalc files out manually would be quite time consuming, so two batch programs were written in Visual C# to automatically write and run numerous ThermoCalc log files.

---

Table 1: Range of 20-32Nb Stainless Steel compositions investigated in the ThermoCalc test matrix. Compositional ranges correspond to ASTM A351, grade CT15C.

	Ni	Cr	Nb	Si	Mn	C	N/Ti	W,Mo,Ti,Zr
Compositional Range (wt%)	31-34	19-21	0.5-1.5	0.5-1.5	0.15-1.5	0.05-0.15	0-0.3	< 0.05

---

The first program titled BatchTCC_Comp shown in Figure 1 takes the compositional range specified, increments each component by the # of iterations set by the user, and then compiles .tcm files for each chemistry combination and runs the scripts. The y-axis is specified by the user where depending on their choice either a stepping or mapping procedure will take place. For this program the x-axis is suppose to be a component, but other properties like temperature and pressure can be used as well. The program outputs both the postscript, and the excel files for the phase or property diagram, and also graphs the composition for each stable phase in the calculation. For this script the relevant phases must be specified, which can be determined by manually doing a test run of the system through ThermoCalc beforehand. The second program titled BatchTCC_Temp is a much simpler version of the first program where only the temperature can be specified as the x-axis for a stepping procedure. If all of the alloy element boxes are used, and the # of iterations is limited to three, the process should be expected to take at least 30 hrs if the ‘# of simultaneous processes’ is set to 3-5. Both the executable installations, and the source code Visual Studio solutions can be downloaded from http://tinyurl.com/9rmjuec.

---

Fig. 1: GUI interface for batch running a compositional matrix through ThermoCalc where the x-axis is a component of the system, and the y-axis is either the phase fraction, another component of the system, or temperature.

---

---

Fig. 2: GUI interface for batch running a compositional matrix through ThermoCalc where the x-axis is a temperature range.

1.1. Preparing for Compilation

After the process has been completed it is important to copy all of the excel file outputs into a separate folder, as some of the subsequent scripts can not handle searching for files through subfolders. Secondly, Thermocalc outputs the excel files in .xls format which conflicts with some of the Matlab code that can only open .xlsx formatted excel files. The first step is to go into the top folder of the directory where the batch ThermoCalc process saved the files. The batch process program saves all of the steping, or mapping plots in the format “< element >< composition >...”, where for a 20Cr32Ni alloy the file name would be “Cr20Ni32...”. For Windows 7 operating systems, find the search bar in the file explorer, and type “Cr*.xls”. This will find all the excel files that start with the characters ‘Cr’. Copy and save all of these files to another folder. Next, download the files named xls2xlsx.rar from http://tinyurl.com/9zsdwf6. Extract, and open the excel file. In the first section type in the path of the copied .xls files in the “Original File Path” cell, and make a new folder for the .xlsx files, and type its path in the “Destination File Path” cell. Click the “Convert to xlsx” button, and all the .xls file should be saved as .xlsx files in the new folder. The same can be done for any compositional data for a specific phase by specifying the cells under the “Compostion xls File” title.

2. Matlab Code - Compiling Mass ThermoCalc Data To Excel

Once all of the ‘.xlx’ files have been converted to ‘.xlsx’ format, and are copied to a single file (eg. ‘_xlsx’), the Matlab distribution files for the batch ThermoCalc application can be employed to extract, and compile all of the relevant information needed to produce matrix arrays, and linear regression fitting. From the volume fraction plots, there are specific points on the curves for each phase that can be useful when trying to quantify the optimization of an alloy based on equilibrium microstructure.

2.1. Finding the solubility temperature, and the terminal phase fractions for intermetallics and Chrome Carbides

The first Matlab script will be used for compiling all the relevant information about any intermetallic phases, or chrome carbides that may have precipitated during long term aging. For example, the G-phase curve in  Figure 3a has three definable points that can be analyzed for all the chemistries when determining the minimization of G-phase. The first definable point is the temperature at which the phase starts to precipitate, referred to as the stability temperature. The chemical formula of G-phase is Ni₁₆Nb₆Si₇, where niobium will be needed to facilitate the precipitation of this intermetallic. The precipitation of G-phase is known to occur through a transformation mechanism with NbC [1, 5, 6]; However, if there is excess niobium at the dendrite boundaries the G-phase will be able to precipitate without consuming NbC precipitates. Excess niobium precipitating out as G-phase is indicative of the plateau region shown in  Figure 3a, where NbC remains unaffected. Once a certain temperature is reached M23C6 becomes stables, promoting the dissolution of NbC, freeing up the niobium to precipitate even more G-phase. Once all of the NbC has been dissolutioned and all the niobium has been exhausted in solution and precipitated as G-phase; the maximum, and terminal phase fraction for G-phase has been reached, denoted as point 3 on  Figure 3a. Point 3 must also be defined as the terminal phase fraction, meaning that for any subsequent decreases in temperature there will be only minor fluctuations in the equilibrium phase fraction. This will avoid any situations where the driving force of the examined phase would start to decrease with decreasing temperature, where it could eventually become metastable, and replaced with another phase.

For Matlab to be able to define these points for each chemistry in the design matrix they must be formulated as mathematical expressions. The stability temperature can be defined as m^α > 0, where m is the phase fraction (mol%) for phase α. Identifying the Plateau region can be achieved by finding where ∂m^α∕∂T ≈ 0, and ∂²m^α∕∂T² ≈ 0. Finding the maximum terminal phase fraction is done in the same way as the plateau region with the added conditions: max(m^α) , and at max(m^α), m^α∕∂T² ≈ 0∀T ∈{T_x…T_rt}, where T_x is the temperature where max(m^α).

(a)

(b)

Fig. 3: The volume fraction of each phase over a temperature range contains three notable regions; the stability temperature, the plateau region, and the maximum volume fraction. These three points are outlined for a) G-Phase, b) and M23C6 of an alloy with the composition Cr19-Ni31-Si0.5-Nb0.5-C0.05-Mn0.825-Ti0.05.

The Matlab distribution files for this section can be found at http://tinyurl.com/8fq4mrj, under ‘precipitate_script’. The Matlab files are given separately, but you can download the full distribution by downloading the ‘ThermoCalc_Script1.rar′ file. The two main scripts that run the program are ‘allxls.m′, and ‘findeq.m′. ‘allxls.m′ is used as a batch processing of all the .xlsx files produced by ThermoCalc, where each file in the directory listed on Line 1 of the script will go through the ‘findeq.m′ script and then appended to a excel file for further use. ‘findeq.m′ is the script that extracts the information described above, shown in  Figure 3a, and  Figure 3b.

The first part of the script extracts the relevant data, and filters out all the unnecessary phases. Line 23 & 24 gets the row index for the row closest to T=400∘C, and T=100∘C. Line 29 extracts the data between this temperature range for the current phase in the for loop. The temperature range 100-400∘C was chosen as the intermetallic phases and the the M23C6 phase have been observed to be stable within this temperature range with ThermoCalc. This is mainly due to the underestimation in the driving force of G-phase predicted by ThermoCalc in the TTNI8 database. Once G-phase is added to an iron database (TTFE*), this range is likely to change to above 800∘C, as G-phase is seen to be stable in niobium stainless steel alloys at this temperature [5, 23, 24]. The If statement in line 31 passes each phase to the next step if the average phase fraction between 100-400K is less than 7%, the maximum phase fraction is less than 10%, and the mean lower threshold is above zero. This is assuming an upper bound on the phases under question, which is a reasonable assumption for the 2032Nb precipitates as literature has never cited a phase fraction for any intermetallics (e.g. G-phase) or chromium carbides above 10 vol%. Setting these bounds also helps filter out strange FCC iterations (eg. FCC_A1#3, FCC_A1#4) that ThermoCalc predicts for very low temperatures, which have a very high driving force. Experimentally no other FCC phases other than MC carbides, and the solution phase were found, so it is appropriate to disregard these other FCC iterations. It should be noted that if this script is to be applied to other systems these conditions may need to be changed. The last conditions filters out the (Nb,Ti)(C,N) phase, or FCC_A1#2 phase, which is addressed in the second, separate script. finphase on line 34 collects the phase names of the phases that pass the criteria in Line 32. The next loop goes through each phase in finphase.

[~,Tx]=min(abs(data(:,1)-400));
[~,Tend]=min(abs(data(:,1)-100));

for phase = phases
 [~,y]=find(ismember(headertext, phase(1,1)));
 %get submatrix of phase from 400C to rt.
 subphase=data(Tx:Tend,y);

 if mean(subphase) < 0.07 && mean(subphase) > 6E-4 ...
	 && max(subphase)< 0.10 ...
	 && strcmp(phase(1,1), 'NP(FCC_A1#2)')==0 
		 finphase = [finphase phase(1,1)];
 end
end

For each phase the derivative, and second derivative of the curve are calculated using the forwards difference method (i → i + 1). The slope between points i, and i + a was also calculated where a is an arbitrary value used to determine if the curve has become constant for further decreases in temperature. For the experiments discussed in this thesis the value of a was set at 8.

% Gets index of current phase for data
[~,y]=find(ismember(headertext, phase(1,1)));
datPhase = data(:,y);
dn = length(datPhase)-1;
dphase = zeros(dn,1);
dphaseI = zeros(dn,1);
ddphase = zeros(dn-1,1);
dChange = 8;

for i=1 : dn-dChange
 dphaseI(i)=abs((datPhase(i+1)-datPhase(i))/(temp(i)-(temp(i+1)-0.1)));
end
for i=1 : dn-dChange
 dphase(i)=abs((datPhase(i+dChange)-datPhase(i))/(temp(i)-temp(i+dChange)));
end
for i=1 : dn-dChange
 ddphase(i)=abs((dphaseI(i+1)-dphaseI(i))/(temp(i)-(temp(i+1)-0.1)));
end
npm = data(:,y);
dtemp = temp(2:length(temp));

Line 63 creates an index n, for counting the number of points in dphase(i), which calculated the derivative of the original curve between i, and i + 8. mxval is the maxium phase fraction attained by the original curve, while mxIdx is the index of this value, and stabTemp is the index of first instance of m^α > 0.

 n = length(dphase)-1;
 k=1;
 mxVal = max(npm);
 mxIdx = find(npm==mxVal, k,'first');
 stabTemp = find(npm,k,'first');

Next, the script needs to go through the phase fraction curve, and its derivative, and find the maximum terminal phase fraction, and if there is a plateau region. This can be achieved by going along each curve, first determining if the driving force for precipitation remains constant, or the phase fraction remains constant, between i and i + a, and then if this point is in a plateau region, or if it is the terminal phase fraction.

if (dphase(k,1) < 2E-5) && (abs(npm(k+1,1)) > 6E-4) ... 
	&& mean(npm(k:k+stbRng,1)) <= mxVal ...
	&& mean(npm(k:k+stbRng,1)) > mxVal-0.003 ...
	&& mean(ddphase(k:k+stbRng,1))< 2E-6...
	&& (dtemp(k) < 800)

	break;
elseif(dphase(k,1) < 2E-5) && (abs(npm(k+1,1)) > 6E-4) ... 
	&& mean(ddphase(k:k+stbRng,1))< 2E-6...
	&& (dtemp(k) < 800) && (dtemp(k) > 250)
	platArrPh = [platArrPh; {data(k-1,y)}];
	platArrT = [platArrT; {dtemp(k)}] ;
	isPlateau = 1; 
end

The If statement in line 78 determines if the driving force remains constant between k (which is just another iterative index similar to i), and k + stbRng, where stbRange is the a parameter, and if the phase fraction is terminal. This will pass if the first derivative, or the slope of the curve is close to zero (2E - 5), and if the phase fraction is within a threshold between mxV al (Maximum phase fraction), and an arbitrary lower threshold, mxV al - 0.003. The other conditions of this If statement determines if the second derivative, or the curvature between k, and k + stbRng (i → i + a) is close to zero, and if the temperature is less than a certain temperature (eg. 800). The T < 800∘C was chosen from examination of ThermoCalc outputs, where G-phase, and M23C6 were never stable at temperatures greater than 800ºC. Once G-Phase is included in the TCFE* database this condition may not be true, as the actual transformation temperature for G-phase in creep resistant stainless steels is 852 ºC[25]. If the phase fraction passes these conditions it will be labeled as ′Stable′, meaning a terminal phase is found, and it will break out of the For loop.

If these conditions are not met the phase fraction will go through another If statement on Line 103, which is almost the same as the If statement on Line 78 without the conditions that the phase fraction needs to be close to the maximum value mxV al. If these conditions pass they will be added to an array platArrPh which will contain the phase fractions for all points on the curve that meet these conditions, and another array platArrT which will contain the temperature. A boolean value called ‘isPlateau′ will be given to this phase. Once the terminal phase fraction is found on Line 78, it will go through a gate labeling the phase as ′PlateauStable′ to indicate that a plateau was found, and the average of platArrPh, and platArrT will be taken and reported as the point where the plateau of the curve is found. If no terminal phase fraction is found, the maximum of the phase will simply be reported, and the phase will be labeled ′NotStable′ if no plateau is found, and ′PlateauNot′ if a plateau is found.

This script will be executed for each composition in the ’allxls.m’ file, and will output the results to a .xlsx file. An example of the resulting table is shown in Table 2. This table can now be reordered by grouping the data based on each phase for comparative purposes.

2.2. Finding Maximum (Nb,Ti)(C,N) Phase fraction, composition, and terminal temperature

The second Matlab script will be used to compile useful information on (Nb,Ti)(C,N) for each chemistry in the compositional matrix. Since (Nb,Ti)(C,N) is of a cubic NaCl structure, which is the same structure as the solution phase (austenite), ThermoCalc creates the (Nb,Ti)(C,N) phase as a new composition set of the FCC phase. Therefore, in ThermoCalc the solution phase is denoted as FCC_A1#1, and (Nb,Ti)(C,N) is denoted as FCC_A1#2. For each chemistry iteration, FCC_A1#2 contains the end members (Nb,Ti)₁(C,N)₁, but the site fractions, or constituents fractions are uncertain. Furthermore, differentiation between NbC, and TiC is not made, and will be grouped into one phase called FCC_A1#2. This information cannot give us distinct phase fractions of NbC, and TiC, as the solubility of titanium in NbC, and niobium in TiC is unknown. For the purpose of this study it will be assumed that these two constituents (niobium, and titanium) are mutually exclusive to their respective phases, and a volume fraction of each will be presented. To determine the composition of the MC carbides, a plot consisting of the mole fraction of constituents of the FCC_A1#2 phase along the temperature range, T_rt → T_m, must be specified.  Figure 4b illustrates an example of a compositional plot of FCC_A1#2, which shows that the composition of FCC_A1#2 can fluctuate based on the solubility of the elements, and the stable phases in the system. To see how the solubilities of titanium and niobium in (Nb,Ti)(C,N) change with variations in alloy composition, the element fractions in FCC_A1#2 can be extracted at specific temperatures for comparison. As illustrated in  Figure 4b Matlab should extract the temperatures and the constituent fraction of FCC_A1#2 at the solubility limit of titanium, the maximum phase fraction of FCC_A1#2, and the dissolution temperature of FCC_A1#2. Due to an error in ThermoCalc reporting constituent fractions of a phase, these fractions do not do not drop to zero once the phase becomes unstable as seen in  Figure 4b.

(a)

(b)

Fig. 4: For the (Nb,Ti)(C) phase, or FCC_A1#2, a) the maximum phase fraction, and the dissolution temperature are extracted, along with b) specific composition of FCC_A1#2 at these regions, and the solubility limit of titanium.

The Matlab distribution files for extracting and compiling the NbC precipitates can be found here at http://tinyurl.com/8fq4mrj, under ‘NbC_script′. Odd anomalies from the composition files makes this script much more particular then the previous script, and should be reviewed before trying to modify it to work with other systems. Some tricks are used to avoid any errors in the data that might be unaccounted for, but some instances will still cause the script to error. If this script is providing trouble, and not completing, it is advised to cut certain sections that might not be needed. For example, there are four parts to this script that extract data from the given .xlsx files. The first part extracts the dissolution temperature and the maximum FCC_A1#2 fraction from the phase fraction data, the second part determines the solubility limit of titanium, the third part takes the maximum phase fraction and determines the respective amounts of NbC and TiC, and the fourth part finds the maximum amount of titanium in FCC_A1#2 and calculates the fractions of NbC and TiC at this point. The fourth part, finding the maximum titanium concentration, is not used in the presented study, but was included for analysis. This fourth part flags most of the errors, so it can be removed if it does not provide critical information to the study.

k=1;
Tx = find(datPhase,k,'first'); % Finds first non-zero
Tend = find(datPhase==0, k,'last'); % Finds last zero for stability range
datPhaseSub = datPhase;
if(datPhase(Tend-1) > 0.01)
 datPhaseSub = datPhaseSub(1:Tend-1);
 Tend = find(datPhaseSub==0, k,'last');
end
if (Tend < Tx+30)
 Tend = length(datPhase);
end

Lines 14-24 determines the suitable range to extract the FCC_A1#2 data. This in an attempt to filter out any noise, or aberrations in the data. T_x finds the first non-zero row, assuming this is the stability temperature for FCC_A1#2, and T_end finds the temperature of the last zero row. Sometimes after dissolution of a phase, ThermoCalc predicts this phase to reprecipitate at low temperatures with a drastically different, often nonsensical compositions. It is important to filter out any regions in the FCC_A1#2 that do not have a composition of (Nb,Ti)(C). Setting T_end to the last zero row will filter out this data, unless FCC_A1#2 drops back down to zero after it has precipitated the second compositionally different phase. The condition on Line 18 will filter out any secondarily precipitated FCC_A1#2 phases by testing if the dissolutioned phase is above 0.01 mole fraction before dissolution, as these new phases usually have a very high driving force. If the FCC_A1#2 is stable until the end of the data, T_end is set to the last temperature in the data as per line 23. Variables dtemp, and npm are extracted between this temperature range (T_x → T_end), which are the temperature array, and the phase fraction array respectively.

mxVal = max(npm);
mxIdx = find(npm==mxVal, k,'first');
mxTemp = dtemp(mxIdx-1);
if(Tend == length(datPhase))
 stabIdx = length(npm);
else
 stabIdx = find(npm==0,k,'first');  
end
stabTemp = dtemp(stabIdx-1);

Lines 31-39 will get the maximum FCC_A1#2 phase fraction, its index, and its corresponding temperature, as well as the dissolution temperature of FCC_A1#2 labelled as stabTemp. Niobium, and titanium mole fraction arrays are extracted and labeled as nbTot, and TiTot. The solubility limit of titanium is found on line 52, where the titanium concentration is below a lower threshold. A While loop introduced on Line 55 is used to find out if the found value is the actual solubility limit of titanium in FCC_A1#2 or if it was picking up some random fluctuation in the data. If the value picked up is a false positive it will be subtracted from the data, and a new solubility limit will be found until the criteria of the While loop is met. The catch statement between lines 59-74 will catch any errors that might occur in this while loop (e.g. indicies out of range), save the error to a log file, and break out of the script, and will continue running with the next file iteration. Lines 31-39 will get the maximum FCC_A1#2 phase fraction, its index, and its corresponding temperature, as well as the dissolution temperature of FCC_A1#2 labelled as stabTemp. Niobium, and titanium mole fraction arrays are extracted and labeled as nbTot, and TiTot. The solubility limit of titanium is found on line 52, where the titanium concentration is below a lower threshold. A While loop introduced on Line 55 is used to find out if the found value is the actual solubility limit of titanium in FCC_A1#2 or if it was picking up some random fluctuation in the data. If the value picked up is a false positive it will be subtracted from the data, and a new solubility limit will be found until the criteria of the While loop is met. The catch statement between lines 59-74 will catch any errors that might occur in this while loop (e.g. indicies out of range), save the error to a log file, and break out of the script, and will continue running with the next file iteration.

% Finds last non-zero for stability range
solIdx = find(tiTot >= 0.5*10^-3, k,'last');
tiTotSub = tiTot;
try
 while(tiTotSub(solIdx-1) < 0.5*10^-3)
  tiTotSub = tiTotSub(1:solIdx-1);
  solIdx = find(tiTotSub >= 0.5*10^-3, k,'last');
 end
catch err
 %open file
 fid = fopen('logFile','a+');
 % write the error to file
 % first line: message
 fprintf(fid,'%s\n',err.message);

 % following lines: stack
 for e=1:length(err.stack)
  fprintf(fid,'%sin %s at %i\n',txt,err.stack(e).name,err.stack(e).line);
 end

 % close file
 fclose(fid)
 return;
end
solIdx = find(tiTot == tiTotSub(solIdx), k,'last');
tempSol = temp(solIdx);   % temp at the sol limit of Ti
phSol = datPhase(solIdx); % phase fraction at the sol limit of Ti.

Lines 79-84 splits the maximum phase fraction found before into NbC, and TiC by finding the concentrations of titanium and niobium in FCC_A1#2, and determining their respective phase fractions, assuming there is no solubility for titanium in NbC, and niobium in TiC.

%Max Phase fraction
tiMx = tiTot(mxIdx);
nbMx = nbTot(mxIdx);
totCon = tiMx + nbMx;
tiPh = tiMx/totCon*mxVal; %TiC Phase fraction
nbPh = nbMx/totCon*mxVal; % NbC Phase fraction

The last part of the script finds the maximum concentration of titanium, checking the validity of the number through a While statement, first checking if the found value is just a large fluctuation in the data, and that the value is below an upper threshold, knowing that the maximum concentration of titanium is never above 0.4. Lines 118-125 determine the respective NbC, and TiC concentrations at the found maximum titanium concentration. An example of the resulting output are shown in Table 3.

%Max Ti
tiTot = tiTot(1:length(tiTot)-1);
mxTi = max(tiTot);
mxTiIdx = find(tiTot==mxTi, k,'first');
tiTotSub = tiTot;
try
 while(abs(mxTi - tiTotSub(mxTiIdx+1)) > 1*10^-3 || mxTi > 0.5)
	if(mxTiIdx < length(tiTotSub)/2)
	 tiTotSub = tiTotSub(mxTiIdx+1:length(tiTotSub));
	else
	 tiTotSub = tiTotSub(1:mxTiIdx-1);
	end
	mxTi = max(tiTotSub);
	mxTiIdx = find(tiTotSub==mxTi, k,'first');
 end
catch err

2.3. Compiling and organizing data in Excel.

Once the Matlab scripts in section 2.1, and section 2.2 have successfully run and output their respective ′.xlsx′ files, download the excel file ′batchTCC.rar′ from this link, http://tinyurl.com/9zsdwf6. Once this file is open, go to the developer tab, and click the Visual Basic Button. In the sortMatlab script, lines 11-12 lines 16-17 refer to the excel files compiled by Matlab. This script will go through the specified files, reorganize, and categorize the data based on the phase and the type of data. The data compiled from section 8.1 is organized in the findPhase function, where the string array strPhase lists the phases that should be searched for. The data array stores the data found for each phase, where this data will be copied to its respective column in the excel table. The phase name, the maximum temperature, and maximum phase fraction, the stability temperature, and any plateau regions, are all stored in the data array. The next function findFCC does the same thing as the findPhase function, but applies it to the second excel file discussed in section 2.2. The organize, and organizeFCC functions take the data arrays returned from findPhase and findFCC, and sort them into their appropriate columns. To run the sortMatlab script press the the sort button on the excel table. Once that has completed the next step is to insert columns after the Name column for each element in the compositional matrix. In the Visual Basic Window click on Model 2 in the Modules folder located on the left-hand side project tree, and run the SplitName script. This takes each string in the Name column and splits the string as to fill the element columns just created. An example of the resulting excel file can be downloaded here, http://tinyurl.com/ccuqjpt. Finally select the header row click on the Data tab in the excel ribbon, and select the Filter option. Each column can now be arranged in ascending or descending order, or with various conditional statements. A supplemental graphing program is included in the .xlsm file, and can be executed by running the openGraphMaker function in Module 2. The X-,Y-,and Z- variables should be selected by clicking the header cells for the columns to be graphed. Up to three constant variables can be selected, by again clicking the header cell of the wanted column, and setting their values to the values you want to keep constant. A text file with the appropriate columns is output for further analysis with other plotting programs (e.g. GLE, or GNUPlot). Figure 4 shows how this plotting tool can be used with the resulting output.