debug & update

debug & update

ReaxFF_Gen/ReaxFF_Gen.m

@@ -11,6 +11,9 @@
 fidin=fopen('input.txt','r');
 fidout=fopen('output.txt','w');
 charset=input('\ncell of elements involved, e.g.{"C","H","O"}:\n');
+example= '{''C'', ''H'', ''O''}';
+fprintf('Please input cell of elements involved, eg.:%s\n',example)
+charset=input('\nPlease input：\n');
 
 while ~feof(fidin)%Number of general parameters
     str=fgetl(fidin);

bonds_analysis/atom_num_autoread.m → bonds_analysis/atom_num_autoreadbond.m

@@ -1,10 +1,10 @@
-%scrit file name atom_num_autoread
+%scrit file name atom_num_autoreadbond
 %purpose:
 %This program is used to read atom number in the simulation box
 %automatically
 %version 1;2024.04.04
 
-function atomnum=atom_num_autoread(dataname)
+function atomnum=atom_num_autoreadbond(dataname)
 fprintf('\nRead atom number from bond.* file, please wait...\n')
 rawdata=fopen(dataname,'r');
 for i=1:3

bonds_analysis/bonds_analysis_speedup.m

@@ -13,7 +13,7 @@
 dataname=input('Please input the file name to be processed: \n','s');
 trajper=input('Please input the output frequency of BO information (Positive integer): \n');
 tartrajectory=input('\nPlease input the timestep of the specified trajectory: \n');
-atomnum=atom_num_autoread(dataname);
+atomnum=atom_num_autoreadbond(dataname);
 disp('bonds_analysis_speedup is running, please wait...')
 tartrajectory={tartrajectory(1)};
 if mod(tartrajectory{1},trajper)~=0
@@ -25,11 +25,9 @@
 end
 
 
-readline=0;
 gap=8+atomnum;
 rawdata=fopen(dataname,'r');
 dataline=fgetl(rawdata);
-readline=readline+1;
 datacell=textscan(dataline,'%s','delimiter','\n');
 datacellchar=char(datacell{1});
 datadel=strrep(datacellchar,'#','');
@@ -39,52 +37,46 @@
     control=0;
 else
     while control
+        gap=8+atomnum;
         i=1;
         unfound=1;
         while unfound
             dataline=fgetl(rawdata);
             while isempty(dataline)
                 dataline=fgetl(rawdata);
             end
-            readline=readline+1;
             i=i+1;
             if i==gap+1
                 unfound=0;
                 break;
             end
         end
-        if mod(readline-1,gap)==0
-            datacell=textscan(dataline,'%s','delimiter','\n');
-            datacellchar=char(datacell{1});
-            datadel=strrep(datacellchar,'#','');
-            datarep=strtrim(datadel);
-            datasplit=strsplit(datarep);
-            if str2num(datasplit{1,2})==tartrajectory{1}
-                control=0;
-            end
-            a0 = ftell(rawdata);
-            dataline=fgetl(rawdata);
-            a1 = ftell(rawdata);
-            dataline=fgetl(rawdata);
-            a2 = ftell(rawdata);
-            datacell=textscan(dataline,'%s','delimiter','\n');
-            datacellchar=char(datacell{1});
-            datadel=strrep(datacellchar,'#','');
-            datarep=strtrim(datadel);
-            datasplit=strsplit(datarep);
-            atomnum = str2double(datasplit{length(datasplit)});
-            fseek(rawdata,-(a2-a0),'cof');
-        else
-            disp('Not timestep row, please check it!!!')
-            return;
+        datacell=textscan(dataline,'%s','delimiter','\n');
+        datacellchar=char(datacell{1});
+        datadel=strrep(datacellchar,'#','');
+        datarep=strtrim(datadel);
+        datasplit=strsplit(datarep);
+        if str2num(datasplit{1,2})==tartrajectory{1}
+            control=0;
         end
+        a0 = ftell(rawdata);
+        dataline=fgetl(rawdata);
+        a1 = ftell(rawdata);
+        dataline=fgetl(rawdata);
+        a2 = ftell(rawdata);
+        datacell=textscan(dataline,'%s','delimiter','\n');
+        datacellchar=char(datacell{1});
+        datadel=strrep(datacellchar,'#','');
+        datarep=strtrim(datadel);
+        datasplit=strsplit(datarep);
+        atomnum = str2double(datasplit{length(datasplit)});
+        fseek(rawdata,-(a2-a0),'cof');
     end
 end
 
 found=6;
 while found
     dataline=fgetl(rawdata);
-    readline=readline+1;
     found=found-1;
 end
 
@@ -98,7 +90,6 @@
 line=2;
 while atomnum
     dataline=fgetl(rawdata);
-    readline=readline+1;
     atomnum=atomnum-1;
     if atomnum<0
         break;

chemi_mechanism/atom_num_autoreadbond.m

@@ -0,0 +1,24 @@
+%scrit file name atom_num_autoreadbond
+%purpose:
+%This program is used to read atom number in the simulation box
+%automatically
+%version 1;2024.04.04
+
+function atomnum=atom_num_autoreadbond(dataname)
+fprintf('\nRead atom number from bond.* file, please wait...\n')
+rawdata=fopen(dataname,'r');
+for i=1:3
+    dataline=fgetl(rawdata);
+end
+datacell=textscan(dataline,'%s','delimiter','\n');
+datacellchar=char(datacell{1});
+datarep=strtrim(datacellchar);
+if contains(datarep,'# Number of particles')
+    fprintf('\nThe atom number line is found, ready to read it\n')
+end
+datadelimiter={'# Number of particles'};
+[C,~]=strsplit(datarep,datadelimiter,'CollapseDelimiters',false);
+atomnum=str2num(C{2});
+fprintf('The atom number is found, which is%d\n',atomnum)
+fclose(rawdata);%
+end

chemi_mechanism/atom_num_autoread.m → chemi_mechanism/atom_num_autoreadtrj.m

@@ -1,10 +1,10 @@
-%scrit file name atom_num_autoread
+%scrit file name atom_num_autoreadtrj
 %purpose:
 %This program is used to read atom number in the simulation box
 %automatically
 %version 1;2023.09.05
 
-function atomnum=atom_num_autoread(datanametrj)
+function atomnum=atom_num_autoreadtrj(datanametrj)
 rawdatatrj=fopen(datanametrj,'r');
 for i=1:3
     dataline=fgetl(rawdatatrj);

chemi_mechanism/bonds_analysis_speedup.m

@@ -21,10 +21,8 @@
 end
 
 rawdata=fopen(datanamebond,'r');
-readline=0;
-gap=8+atomnum;
+atomnum=atom_num_autoreadbond(datanamebond);
 dataline=fgetl(rawdata);
-readline=readline+1;
 datacell=textscan(dataline,'%s','delimiter','\n');
 datacellchar=char(datacell{1});
 datadel=strrep(datacellchar,'#','');
@@ -34,58 +32,52 @@
     control=0;
 else
     while control 
+        gap=8+atomnum;
         i=1;
         unfound=1;
         while unfound
             dataline=fgetl(rawdata);
             while isempty(dataline)%可能存在ITEM: TIMESTEP行前空行
                 dataline=fgetl(rawdata);
             end
-            readline=readline+1;
             i=i+1;
             if i==gap+1
                 unfound=0;
                 break;
             end
         end
-        if mod(readline-1,gap)==0
-            if feof(rawdata)
-                fprintf('\nThe last line of bonds.* file is reached. Please check if the atom number or timestep is correct.\n The timestep in species.* and bond.* maybe inconsistent, eg. 29,429...\n')
-                error('As above, please check it!')
-            end
-            datacell=textscan(dataline,'%s','delimiter','\n');
-            datacellchar=char(datacell{1});
-            datadel=strrep(datacellchar,'#','');
-            datarep=strtrim(datadel);
-            datasplit=strsplit(datarep);
-            if str2num(datasplit{1,2})==tartrajectory{1}
-                control=0;
-            elseif str2num(datasplit{1,2})>tartrajectory{1}
-                error('The frame No. appointed according to species.* file exceeds that in bonds.* file, please check it!')
-            end
-            a0 = ftell(rawdata);
-            dataline=fgetl(rawdata);
-            a1 = ftell(rawdata);
-            dataline=fgetl(rawdata);
-            a2 = ftell(rawdata);
-            datacell=textscan(dataline,'%s','delimiter','\n');
-            datacellchar=char(datacell{1});
-            datadel=strrep(datacellchar,'#','');
-            datarep=strtrim(datadel);
-            datasplit=strsplit(datarep);
-            atomnum = str2double(datasplit{length(datasplit)});
-            fseek(rawdata,-(a2-a0),'cof');
-        else
-            disp('This is not a line with timestep, please check it!!!')
-            return;
+        if feof(rawdata)
+            fprintf('\nThe last line of bonds.* file is reached. Please check if the atom number or timestep is correct.\n The timestep in species.* and bond.* maybe inconsistent, eg. 29,429...\n')
+            error('As above, please check it!')
+        end
+        datacell=textscan(dataline,'%s','delimiter','\n');
+        datacellchar=char(datacell{1});
+        datadel=strrep(datacellchar,'#','');
+        datarep=strtrim(datadel);
+        datasplit=strsplit(datarep);
+        if str2num(datasplit{1,2})==tartrajectory{1}
+            control=0;
+        elseif str2num(datasplit{1,2})>tartrajectory{1}
+            error('The frame No. appointed according to species.* file exceeds that in bonds.* file, please check it!')
         end
+        a0 = ftell(rawdata);
+        dataline=fgetl(rawdata);
+        a1 = ftell(rawdata);
+        dataline=fgetl(rawdata);
+        a2 = ftell(rawdata);
+        datacell=textscan(dataline,'%s','delimiter','\n');
+        datacellchar=char(datacell{1});
+        datadel=strrep(datacellchar,'#','');
+        datarep=strtrim(datadel);
+        datasplit=strsplit(datarep);
+        atomnum = str2double(datasplit{length(datasplit)});
+        fseek(rawdata,-(a2-a0),'cof');
     end
 end
 
 found=6;
 while found
     dataline=fgetl(rawdata);
-    readline=readline+1;
     found=found-1;
 end
 
@@ -99,7 +91,6 @@
 line=2;atomnumcopy=atomnum;
 while atomnumcopy
     dataline=fgetl(rawdata);
-    readline=readline+1;
     atomnumcopy=atomnumcopy-1;
     if atomnumcopy<=0
         break;
@@ -146,6 +137,6 @@
     bondoutdata(:,16:end)=[];
 end
 
-fprintf('\nbonds_analysis_speedup is successfully finished, BO information is saved in bondoutdata, search line number is recorded in readline\n')
+fprintf('\nbonds_analysis_speedup is successfully finished, BO information is saved in bondoutdata\n')
 clear atomnumcopy ans bondnumdata control datacell datacellchar datadel dataline datarep datasplit found gap i j k kk line 
 clear outputans unfound dataoutrow dataoutcol dataoutputrow dataoutcolchar dataoutputcol filename a0 a1 a2

chemi_mechanism/chemi_mechanism.m

@@ -29,7 +29,7 @@
 trajper=input('\neg., [400 400 400]:  \n');
 
 fprintf('Automatically read atom number from the *.lammpstrj file, please wait...')
-atomnum=atom_num_autoread(datanametrj); 
+atomnum=atom_num_autoreadtrj(datanametrj); 
 
 check_control_origin=input('\nPlease input the number for frame No. check of bonds.* and *.lammpstrj files compared with species file，\navoiding mismatch induced error, >=5 is suggested: \n');
 if check_control_origin <4

chemi_mechanism/lammpstrj_analysis.m

@@ -18,13 +18,11 @@
 end
 fprintf('\nlammpstrj_analysis is running, please wait...')
 
-readline=0;
+atomnum=atom_num_autoreadtrj(datanametrj);
 gap=9+atomnum;
 
 dataline=fgetl(rawdatatrj);
-readline=readline+1;
 dataline=fgetl(rawdatatrj);
-readline=readline+1;
 datacell=textscan(dataline,'%s','delimiter','\n');
 datacellchar=char(datacell{1});
 datarep=strtrim(datacellchar);
@@ -40,49 +38,44 @@
     return;
 elseif str2num(datarep)<tartrajectory{1}
     while control
+        gap=9+atomnum;
         i=1;
         unfound=1;
         while unfound
             dataline=fgetl(rawdatatrj);
             while isempty(dataline)%可能存在ITEM: TIMESTEP行前空行
                 dataline=fgetl(rawdatatrj);
             end
-            readline=readline+1;
             i=i+1;
             if i==gap+1
                 unfound=0;
                 break;
             end
         end
-        if mod(readline-2,gap)==0
-            if feof(rawdatatrj)
-                fprintf('\nThe last line of lammpstrj.* file is reached. Please check if the atom number or timestep is correct.\n The timestep in species.* and lammpstrj.* maybe inconsistent, eg. 29,429...\n')
-                error('As above, please check it!')
-            end
-            datacell=textscan(dataline,'%s','delimiter','\n');
-            datacellchar=char(datacell{1});
-            datarep=strtrim(datacellchar);
-            if str2num(datarep)==tartrajectory{1}
-                control=0;
-            elseif str2num(datarep)>tartrajectory{1}
-                error('The frame No. appointed according to species.* file exceeds that in lammpstrj.* file, please check it!')
-            end
-            a0 = ftell(rawdatatrj);
-            dataline=fgetl(rawdatatrj);
-            a1 = ftell(rawdatatrj);
-            dataline=fgetl(rawdatatrj);
-            a2 = ftell(rawdatatrj);
-            datacell=textscan(dataline,'%s','delimiter','\n');
-            datacellchar=char(datacell{1});
-            datadel=strrep(datacellchar,'#','');
-            datarep=strtrim(datadel);
-            datasplit=strsplit(datarep);
-            atomnum = str2double(datasplit{length(datasplit)});
-            fseek(rawdatatrj,-(a2-a0),'cof');
-        else
-            disp('This is not a line with timestep information, please check it!!!')
-            return;
+        if feof(rawdatatrj)
+            fprintf('\nThe last line of lammpstrj.* file is reached. Please check if the atom number or timestep is correct.\n The timestep in species.* and lammpstrj.* maybe inconsistent, eg. 29,429...\n')
+            error('As above, please check it!')
+        end
+        datacell=textscan(dataline,'%s','delimiter','\n');
+        datacellchar=char(datacell{1});
+        datarep=strtrim(datacellchar);
+        if str2num(datarep)==tartrajectory{1}
+            control=0;
+        elseif str2num(datarep)>tartrajectory{1}
+            error('The frame No. appointed according to species.* file exceeds that in lammpstrj.* file, please check it!')
         end
+        a0 = ftell(rawdatatrj);
+        dataline=fgetl(rawdatatrj);
+        a1 = ftell(rawdatatrj);
+        dataline=fgetl(rawdatatrj);
+        a2 = ftell(rawdatatrj);
+        datacell=textscan(dataline,'%s','delimiter','\n');
+        datacellchar=char(datacell{1});
+        datadel=strrep(datacellchar,'#','');
+        datarep=strtrim(datadel);
+        datasplit=strsplit(datarep);
+        atomnum = str2double(datasplit{length(datasplit)});
+        fseek(rawdatatrj,-(a2-a0),'cof');
     end
 end
 %
@@ -172,7 +165,6 @@
 trjdata=[];line=1;atomnumcopy=atomnum;
 while atomnumcopy
     dataline=fgetl(rawdatatrj);
-    readline=readline+1;
     atomnumcopy=atomnumcopy-1;
     if atomnumcopy<=0
         break;