ICME 2012 POT-GEN-3
Author(s): Laalitha Liyanage
Contents |
License
By using the codes provided here you accept the the Mississippi State University's license agreement. Please read the agreement carefully before usage.
Overview
Potential space sampling script will sample the parameter space specified by user (sensitive parameters and bounds) utilizing a stratified sampling method known as Latin Hypercube Sampling.
Modifying the script
- The parameters which are sensitive and their ranges should be specified at the beginning of the script. As an example all parameters are specified with the ranges used in script-2.
- The number of samples can be changed to a desired value. A small value will represent the parameter space less accurately compared to a larger value. Also a larger value means more calculations and more time requirements.
Running the script
- Save the script to LINUX text file and name it any_name.m.
- Copy and save auxiliary scripts with provided names.
The script can be executed in a terminal or using a pbs command script using the following command
matlab -nodisplay -nosplash < any_name.m
The script will call script-1. Therefore the number of processors allocated in the pbs command script should correspond with the number of processors specified in the mpirun command in script-1.
Potential Space Sampling script
%% Sampling the potential space
clear all
%Specify sensitive parameters and their bounds
%Example specifies all parameters and maximum ranges
%X(1:12) = [alpha, asub, b0 , b1 , b2 , b3 , t0 , t1 , t2 ,t3 , cmax, cmin]
low = [3.78 , 0.8, 0.0, 0.0, 0.0, 0.0 , -10, -10 , -10, -10, 2.0, 0.0];
high = [5.67, 1.2, 10, 10, 10, 10 , 10, 10 , 110, 10, 2.8, 2.0];
% Specify no. of samples
ndiv = 2000;
% Specify no. of sensitive parameters
npar = 12;
% Design LHS
X = lhsdesign(ndiv,npar,'criterion','correlation');
%Initialize database for responses
Y = zeros([size(X,1),15]);
%% Initialize potential parameters to default in 01_Al_Coarse_pot_eval.m
%Expermental data
Ec = 3.43;
%DFT data - Data from your calculations of E-V curves
a0 = 4.05;
B = 82/160.21765;
V0 = a0^3/4;
%calculate alpha
alpha = sqrt(9*B*V0/Ec);
%X(1:12) = [alpha, asub, b0 , b1 , b2 , b3 , t0 , t1 , t2 , t3 , cmax, cmin]
default = [alpha, 1.07, 2.21, 2.20, 6.0, 2.2 , 1.0, -1.78, -2.21, 8.01, 2.8, 2.0];
clear parameter
parameter.elt = 'Al';
parameter.lat = 'fcc';
parameter.z = 12.0;
parameter.ielement = 13.0;
parameter.atwt = 26.9815;
parameter.rozero = 1.0;
parameter.ibar = 0;
parameter.rcut = 5.0;
parameter.alat = a0;
parameter.esub = Ec;
parameter.alpha = default(1);
parameter.asub = default(2);
parameter.b0 = default(3);
parameter.b1 = default(4);
parameter.b2 = default(5);
parameter.b3 = default(6);
parameter.t0 = default(7);
parameter.t1 = default(8);
parameter.t2 = default(9);
parameter.t3 = default(10);
parameter.cmax = default(11);
parameter.cmin = default(12);
%load('parameter.mat')
%modify parameter
filename1 = 'library.meam.init';
filename2 = 'Al.meam.init';
write_lmpmeam(filename1,filename2,parameter);
%% Run LAMMPS and get responses for sample potentials
filename1 = 'library.meam';
filename2 = 'Al.meam';
for i = 1:size(X,1)
p = parameter;
% Make changes
p.alpha = low(1)+(high(1)-low(1))*X(i,1);
p.asub = low(1)+(high(2)-low(2))*X(i,2);
p.b0 = low(1)+(high(3)-low(3))*X(i,3);
p.b1 = low(1)+(high(4)-low(4))*X(i,4);
p.b2 = low(1)+(high(5)-low(5))*X(i,5);
p.b3 = low(1)+(high(6)-low(6))*X(i,6);
p.t0 = low(1)+(high(7)-low(7))*X(i,7);
p.t1 = low(1)+(high(8)-low(8))*X(i,8);
p.t2 = low(1)+(high(9)-low(9))*X(i,9);
p.t3 = low(1)+(high(10)-low(10))*X(i,10);
p.cmax = low(1)+(high(11)-low(11))*X(i,11);
p.cmin = low(1)+(high(12)-low(12))*X(i,12);
write_lmpmeam(filename1,filename2,p);
system('./bb_script.py');
p_Ec = find_energy('responses.txt', 'E_coh');
p_a0 = find_energy('responses.txt', 'a_0');
p_V0 = find_energy('responses.txt', 'V_0');
p_c11 = find_energy('responses.txt', 'C11');
p_c12 = find_energy('responses.txt', 'C12');
p_c44 = find_energy('responses.txt', 'C44');
p_B = find_energy('responses.txt', 'B');
p_E_100 = find_energy('responses.txt', 'E_100');
p_E_110 = find_energy('responses.txt', 'E_110');
p_E_111 = find_energy('responses.txt', 'E_111');
p_GSFE_I = find_energy('responses.txt', 'GSFE_I');
p_GSFE_E = find_energy('responses.txt', 'GSFE_E');
p_Evac = find_energy('responses.txt', 'E_vac');
p_Eoct = find_energy('responses.txt', 'E_Oct');
p_Etet = find_energy('responses.txt', 'E_Tetra');
Y(i,:) = [p_Ec, p_a0, p_V0, p_c11, p_c12, p_c44, p_B, p_E_100, p_E_110,...
p_E_111, p_GSFE_I, p_GSFE_E, p_Evac, p_Eoct, p_Etet];
end
save('Al_LHS_v1.mat','X','Y','default','low','high')
%{
%% Further analysis...
clear all
load('Al_LHS_v1.mat')
% Select potential based on criteria and objective function
%Set goal values at the beginning of Potential Space Evaluation Script
load('goal.mat'); %Load values for goal
%Example weight
w = [1,1,1,1,1,1,1,1,1,1,1,1,1,1,1];
clear fobj
for i = 1:size(Y,1)
fobj(i,:) = w.*((Y(i),:) - goal).^2;
end
fobj = sum(fobj,2);
nmin = find(fobj == min(fobj));
X(nmin(1),:);
Y(nmin(1),:);
% Make New potential files from X(i,:)
i=nmin(1);
load('parameter.mat')
p = parameter;
p.alpha = low(1)+(hi(1)-low(1))*X(i,1);
p.asub = low(1)+(hi(2)-low(2))*X(i,2);
p.b0 = low(1)+(hi(3)-low(3))*X(i,3);
p.b1 = low(1)+(hi(4)-low(4))*X(i,4);
p.b2 = low(1)+(hi(5)-low(5))*X(i,5);
p.b3 = low(1)+(hi(6)-low(6))*X(i,6);
p.t0 = low(1)+(hi(7)-low(7))*X(i,7);
p.t1 = low(1)+(hi(8)-low(8))*X(i,8);
p.t2 = low(1)+(hi(9)-low(9))*X(i,9);
p.t3 = low(1)+(hi(10)-low(10))*X(i,10);
p.cmax = low(1)+(hi(11)-low(11))*X(i,11);
p.cmin = low(1)+(hi(12)-low(12))*X(i,12);
write_lmpmeam('library.meam.new','Al.meam.new',p);
%}
OUTPUT of script
The script will output a data file Al_LHS_v1.mat that will contain the samples and their respective responses and the bounds of the parameters used. This file is used by the script-4 to generate analytical models of the response surface.
Further Analysis
The further analysis section of the script includes procedures to do a multiobjective optimization on the samples generated and look at an initial optimized potential. This section can be run via the MATLAB GUI. Do swsetup matlab and type matlab. Open script highlight everthing between %{ and '%} after the save command. Then hit F9 key. A set of values for the parameters and its responses should come up.
