Welcome to Liren Liu’s site!¶
About Me:¶
Liren Liu (刘立仁)¶
Research Interests: Mechanisms of 2D materials growth;Machine learning method;Condensed matter science;Cluster Science.
Research Skills: VASP;GAUSSIAN;LAMMPS;DFTB+;PYTHON;LIUNX
Education and work experience:
2017.7– ,Nanjing University of Aeronautics and Astronautics, Postdoctal (Zhuhua Zhang)
2013.9-2017.7, University of Science and Technology of China,Doctor (Jinlong Yang)
Publications¶
1.Yongmin He,Liren Liu#,Chao Zhu,Prafful Golani,Pengyi Tang, Caitian Gao, Iordi Arbiol,Qijie Wang,Zhuhua Zhang*,Zheng Liu*. Reproducibility of atom-thin catalyst in layered platinum sulfide,Science,submitting.2019
2.Liren Liu,Zhuhua Zhang*,Xiolong Liu,Xiaoyu Xuan,Mark C.Hersam*,Boris I.Yakobson and Wanlin Guo*.Borophene concentric superlattices via self-assembly of twin boundaries.Nat. Comm., submitting.2019
3.Longjuan Kong,Liren Liu#,Lan Chen,Qing Zhong,Peng Cheng,Hui Li,Zhuhua Zhang*,Kehui Wu*.One-dimensional nearly free electron states in borophene,Nanoscal,11,15605.2019
4.Zhili Hu, Zhuhua Zhang, Liren Liu, Wanlin Guo. Extreme pseudomagnetic fields in carbon nanocones by simple loads. J. Mech. Phys. Solids,12,1-9.2019
5.Lijuan Yan,Jiamei Shao,Liren Liu*,Chunlei Chen*.Seventeen-coordinate binary metal superatoms: M@Li17,Chem.Phys.Chem.Lett,733,136693.2019
6.Miaomiao Ma,Liren Liu*,Hengjiang Zhu*,Junzhe Lu,Guipin Tan.Structural evolution and properties of small-size thiol-protected gold nanoclusters.Mole.Phys.116,1804.2018
7.Liren Liu,Longjiu Cheng,Jinlong Yang.The superatom molecule theory of metal clusters.Sci.Sin. Chim.48,00167.2018
8.Jishi Chen, Liren Liu#, Xu Liu, Lingwen Liao, Shengli Zhuang, Shiming Zhou, Jinlong Yang*, Zhikun Wu*.Gold-doping of Double-crown Pd Nanoclusters. Chem.Eur.J.,23,18187–18192.2017
9.Liren Liu, Jinyun Yuan, Longjiu Cheng, Jinlong Yang*. New insights into the stability and structural evolution of some gold nanoclusters. Nanoscale,9, 856-861.2017
10.Liren Liu, Pai Li, Lan-Feng Yuan*, Longjiu Cheng, Jinlong Yang*. From Isosuperatoms to Isosupermolecules: New Concepts in Cluster Science. Nanoscale, 8, 12787-12792.2016
11.Liren Liu,Yanbo Zou and Hengjiang Zhu*.Structure and electronic properties of GaN tubelike clusters and single-walled GaN nanotubes.Int.J.Mode.Phys.B.,29,1550116.2015
12.Jishi Chen, Liren Liu#, Linhong Weng, Yuejian Lin, Lingwen Liao, Chengming Wang, Jinlong Yang*, Zhikun Wu*. Synthesis and properties evolution of a family of tiara-like phenylethanethiolated palladium nanoclusters.Sci.Rep.,5:16628.2015
13.Chuanhao Yao,Jishi Chen,Man-Bo Li,Liren Liu,Jinlong Yang,Zhikun Wu*. Adding two active silver atoms on Au25 nanoparticle.Nano Lett., 15 (2), pp 1281–1287.2015
14.Lingwen Liao,Shiming Zhou,Yafei Dai,Liren Liu,Chuanhao Yao,Cenfeng Fu,Jinlong Yang*,Zhikun Wu*. Mono-Mercury Doping of Au25 and the HOMO/LUMO Energies Evaluation Employing Differential Pulse Voltammetry. J. Am. Chem. Soc., 137 (30), pp 9511–9514.2015
15.Liren Liu,Xueling Lei,Hang Chen,Hengjiang Zhu*.Geometry and electronic properties of Bn(n=2-15) clusters.Chin.Phys.Soc.,58,5355.2012
Note: Co-first auther(#);Corresponding author(*)
Scripts and codes:¶
Vasp¶
Vasp is a a computer program for atomic scale materials modelling, e.g. electronic structure calculations and quantum-mechanical molecular dynamics, from first principles.
This page contains scripts for vasp post treatment.
Useful tools on web¶
- VTST (Transition State Tools for VASP)
- It is recommended to refer to lipai’s blog for NEB calculation.
- VTSTscripts
- VASPKIT
- p4vasp
- Pymatgen (Python Materials Genomics)
Scripts¶
cif2pos.sh¶
This page contains a bash script to convert cif format of structure into poscar format
#!/bin/sh
input=$1
format=$2
if test -z $2; then
format="ms"
fi
awk 'BEGIN{OFS="\t\t"}{if(NR==1) {print $0;print "1"}
if($1=="_cell_length_a") print $2,"0","0";
if($1=="_cell_length_b") print "0",$2,"0";
if($1=="_cell_length_c") print "0","0",$2;}' $input >POSCAR-cif
awk 'BEGIN{ORS=""}/_atom_site_occupancy/,/loop_/{
if($2~/^[A-Z]+/) {
if(a[i]!=$2) a[++i]=$2
++b[i]}
}END{
for(j=1;j<=i;j++) print a[j],"\t";
print "\n"
for(j=1;j<=i;j++) print b[j],"\t";
print "\n"
}' $input >>POSCAR-cif
echo "Direct" >>POSCAR-cif
awk -v format="$format" 'BEGIN{OFS="\t\t"}/_atom_site_fract_z/,/loop_/{
if($2~/^[A-Z]+/) {
if(format=="vtst") print $4,$5,$6
else if(format=="ms") print $3,$4,$5}
}' $input >>POSCAR-cif
converg.sh¶
This page contains a bash script to supervise the convergence of force in VASP calculation. Fixed atoms are not considered here.
#!/bin/sh
# eliminate forces of fixed atoms
#lipai@mail.ustc.edu.cn
fix=$1
if test -z $1; then
fix=0
fi
# get data form OUTCAR
awk -v fix="$fix" '/POSITION/,/drift/{
if($1~/^[0-9.]+$/&&$3>=fix) print $1,$2,$3,sqrt($4*$4+$5*$5+$6*$6i);
else if($1=="total") print $0
}' OUTCAR >temp.f
awk '{
if($1=="total") {print ++i,a;a=0}
else {if(a<$4) a=$4}
}' temp.f >force.conv
sed -i '1,9d' force.conv
rm temp.f
awk '/entropy=/{if(i<9) i++;else print ++i,$7}' OUTCAR >energy.conv
# plot
gnuplot <<EOF
set grid
set term post
set output 'b.ps'
set xlabel 'Ion steps'
set title 'Energy & Max Force of each ion steps'
set key box reverse spacing 3.0
set ytics nomirror
set y2tics
set ylabel 'Energy(eV)'
set y2label 'Max Force(eV/Angstrom)'
plot 'energy.conv' u 1:2 w l lw 3 lc rgb "red" axes x1y1 t "Energy ",'force.conv' u 1:2 w l lw 2 lc 3 axes x1y2 t "Max Force "
EOF
gs -sDEVICE=jpeg -r300 -sPAPERSIZE=a4 -dBATCH -dNOPAUSE -sOutputFile=conv.jpg b.ps
convert -rotate 90 conv.jpg conv.jpg
mogrify -trim conv.jpg
gnuplot <<EOF
set term dumb
plot 'energy.conv' w l t " Energy "
plot 'force.conv' w l t " Force "
EOF
rm b.ps # force.conv energy.conv
get_vibration.sh¶
This bash script extract vibration info from OUTCAR, which can be further used by jmol for visualization.
You may also refer to my blog for jmol and vibration visualization
#!/bin/sh
#edge by lipai@mail.ustc.edu.cn
pos2xyz.pl POSCAR
n=`grep PiTHz OUTCAR |wc -l`
awk '/PiTHz/,/POTIM/{print $4" "$5" "$6}' OUTCAR> tmp
awk 'BEGIN{RS="\n \n"}{a++}{print >"tmp_"a}' tmp
for i in `seq 1 $n`
do
sed -i "1,2s/.*/ /" tmp_$i
paste -d " " POSCAR.xyz tmp_$i >vib_$i.xyz
done
rm tmp*
mkdir vib
mv vib_* vib
plot_structure.sh¶
This bash script plot POSCAR and CONTCAR in folders ini and fin. If you have many structures to plot, this script can be efficient after necessary modification.
# get picture of ini&fin(POSCAR&CONTCAR)
for i in ini fin
do for j in POSCAR CONTCAR
do pos2xyz.pl $i/$j
cat >xcrysden.script <<EOF
::scripting::open --xyz $i/$j.xyz
# ------------------------------------------------------------------------
# END: XSF structure data
# ------------------------------------------------------------------------
# ======================================================================== #
# #
# STATE-PART OF THE FILE #
# #
# ======================================================================== #
# ------------------------------------------------------------------------
# definition of xcMisc array
# ------------------------------------------------------------------------
array set xcMisc {ImageMagick.import /usr/bin/import reduce_to {} rescale_image_list {up down left right center rotXY rotXZ rotYZ rotAB rotAC rotBC wireframes_2d pointlines_2d pipeballs_2d ballsticks2_2d ballsticks1_2d spacefills_2d spacefills_3d ballsticks_3d pipeballs_3d sticks_3d dm_wire dm_solid dm_anaglyph dm_stereo dm_smooth dm_flat rep_unit rep_asym} resolution 1440x870 gif_encoder convert status_init_label {Building GUI ...} wm_rootYshift 0 movie_encoder ppmtompeg debug 0 ppmtompeg /usr/bin/ppmtompeg convert /usr/bin/convert titlefile /tmp/xc_20968/a.xsf xwd /usr/bin/xwd resolution_ratio1 1.0 resolutionX 1440 resolution_ratio2 1.0 wm_rootXshift 0 resolutionY 870 ImageMagick.convert /usr/bin/convert}
# ------------------------------------------------------------------------
# display "waiting" toplevel and watch-cursor
# ------------------------------------------------------------------------
set wait_window [DisplayUpdateWidget "Reconstructing" "Reconstructing the structure and display parameters. Please wait"]
SetWatchCursor
set xcCursor(dont_update) 1
# ------------------------------------------------------------------------
# size of the main window fonts
# ------------------------------------------------------------------------
wm geometry . 1440x851
# ------------------------------------------------------------------------
# BEGIN: create fonts
# ------------------------------------------------------------------------
saveState:fontCreate font10 -family helvetica -size 10 -weight normal -slant roman -underline 0 -overstrike 0
saveState:fontCreate font1 -family helvetica -size 10 -weight normal -slant roman -underline 0 -overstrike 0
saveState:fontCreate font11 -family helvetica -size 10 -weight normal -slant roman -underline 0 -overstrike 0
saveState:fontCreate font2 -family helvetica -size 8 -weight normal -slant roman -underline 0 -overstrike 0
saveState:fontCreate font12 -family helvetica -size 10 -weight normal -slant roman -underline 0 -overstrike 0
saveState:fontCreate font3 -family helvetica -size 12 -weight normal -slant roman -underline 0 -overstrike 0
saveState:fontCreate font13 -family helvetica -size 10 -weight normal -slant roman -underline 0 -overstrike 0
saveState:fontCreate font4 -family helvetica -size 10 -weight normal -slant roman -underline 0 -overstrike 0
saveState:fontCreate TkDefaultFont -family helvetica -size 10 -weight normal -slant roman -underline 0 -overstrike 0
saveState:fontCreate TkMenuFont -family helvetica -size 10 -weight normal -slant roman -underline 0 -overstrike 0
saveState:fontCreate font14 -family courier -size 10 -weight normal -slant roman -underline 0 -overstrike 0
saveState:fontCreate font5 -family helvetica -size 10 -weight normal -slant roman -underline 0 -overstrike 0
saveState:fontCreate font6 -family helvetica -size 8 -weight normal -slant roman -underline 0 -overstrike 0
saveState:fontCreate font7 -family helvetica -size 10 -weight normal -slant roman -underline 0 -overstrike 0
saveState:fontCreate TkHeadingFont -family helvetica -size 10 -weight bold -slant roman -underline 0 -overstrike 0
saveState:fontCreate font8 -family helvetica -size 10 -weight normal -slant roman -underline 0 -overstrike 0
saveState:fontCreate font9 -family helvetica -size 10 -weight normal -slant roman -underline 0 -overstrike 0
saveState:fontCreate TkTooltipFont -family helvetica -size 8 -weight normal -slant roman -underline 0 -overstrike 0
saveState:fontCreate TkTextFont -family helvetica -size 10 -weight normal -slant roman -underline 0 -overstrike 0
saveState:fontCreate TkCaptionFont -family helvetica -size 11 -weight bold -slant roman -underline 0 -overstrike 0
saveState:fontCreate TkSmallCaptionFont -family helvetica -size 8 -weight normal -slant roman -underline 0 -overstrike 0
saveState:fontCreate TkFixedFont -family courier -size 10 -weight normal -slant roman -underline 0 -overstrike 0
saveState:fontCreate TkIconFont -family helvetica -size 10 -weight normal -slant roman -underline 0 -overstrike 0
# ------------------------------------------------------------------------
# END: create fonts
# ------------------------------------------------------------------------
# ------------------------------------------------------------------------
# BEGIN: take care of display-mode
# ------------------------------------------------------------------------
set translationStep 0.05
set rotstep 10
set light On
Lighting On
array set mode2D {PL Off SF Off WF Off BS1 Off PB Off BS2 Off}
array set mode3D {space Off sticks On pipe Off balls On}
array set dispmode {mode3D_name BS mode3D Preset mode3D_f2_packinfo {-in .ctrl.c.f.fr3 -anchor center -expand 1 -fill x -ipadx 0 -ipady 0 -padx 0 -pady 0 -side top} style 3D mode3D_ModeFrame .ctrl.c.f.fr3.2.a0}
saveState:displayMode
set style3D(draw) solid; Style3D draw solid
set style3D(shade) smooth; Style3D shade smooth
set viewer(rot_zoom_button_mode) Discrete; Viewer:rotZoomButtonMode
# ------------------------------------------------------------------------
# END: take care of display-mode
# ------------------------------------------------------------------------
# ------------------------------------------------------------------------
# Number of Units Drawn
# ------------------------------------------------------------------------
set nxdir 1
set nydir 1
set nzdir 1
# ------------------------------------------------------------------------
# BEGIN: Atomic-Labels/Fonts
# ------------------------------------------------------------------------
array set atomLabel {atomFont {} globalFont {} atomFont.brightColor {1.0 1.0 1.0} globalFont.brightColor {1.0 1.0 1.0} atomFont.do_display 1 atomFont.label {} globalFont.do_display 1 atomFont.darkColor {0.0 0.0 0.0} atomFont.id {} globalFont.darkColor {0.0 0.0 0.0} fontBrowser {Simple Font Browser}}
set t [ModAtomLabels]
.mesa xc_setfont {} {1.0 1.0 1.0} {0.0 0.0 0.0}
ModAtomLabels:advancedCheckButton default
ModAtomLabels:advancedCheckButton custom
ModAtomLabels:advancedCloseUpdate dummy update
CancelProc \$t
# ------------------------------------------------------------------------
# END: Atomic-Labels/Fonts
# ------------------------------------------------------------------------
# ------------------------------------------------------------------------
# BEGIN: Various colors
# ------------------------------------------------------------------------
xc_newvalue .mesa 8 0 0.770000 1.000000 0.420000
xc_newvalue .mesa 8 1 0.000000 0.950000 0.950000
xc_newvalue .mesa 8 2 0.950000 0.950000 0.950000
xc_newvalue .mesa 8 3 0.950000 0.950000 0.950000
xc_newvalue .mesa 8 4 0.707900 0.707900 0.707900
xc_newvalue .mesa 8 5 0.707900 0.707900 0.707900
xc_newvalue .mesa 8 6 0.950000 0.950000 0.000000
xc_newvalue .mesa 8 7 0.644500 0.804700 0.856900
xc_newvalue .mesa 8 8 0.700000 0.000000 0.000000
xc_newvalue .mesa 8 9 0.834500 0.950000 0.950000
xc_newvalue .mesa 8 10 0.950000 0.950000 0.950000
xc_newvalue .mesa 8 11 0.000000 0.950000 0.950000
xc_newvalue .mesa 8 12 0.000000 0.950000 0.950000
xc_newvalue .mesa 8 13 0.950000 0.000000 0.950000
xc_newvalue .mesa 8 14 0.000000 0.950000 0.950000
xc_newvalue .mesa 8 15 0.950000 0.950000 0.000000
xc_newvalue .mesa 8 16 0.950000 0.950000 0.450000
xc_newvalue .mesa 8 17 0.710000 1.000000 0.000000
xc_newvalue .mesa 8 18 0.950000 0.950000 0.950000
xc_newvalue .mesa 8 19 0.000000 0.950000 0.950000
xc_newvalue .mesa 8 20 0.000000 0.950000 0.950000
xc_newvalue .mesa 8 21 0.000000 0.950000 0.950000
xc_newvalue .mesa 8 22 0.600000 0.600000 0.600000
xc_newvalue .mesa 8 23 0.644500 0.804700 0.856900
xc_newvalue .mesa 8 24 0.644500 0.804700 0.856900
xc_newvalue .mesa 8 25 0.644500 0.804700 0.856900
xc_newvalue .mesa 8 26 0.950000 0.000000 0.000000
xc_newvalue .mesa 8 27 0.644500 0.804700 0.856900
xc_newvalue .mesa 8 28 0.644500 0.804700 0.856900
xc_newvalue .mesa 8 29 0.820000 0.450000 0.140000
xc_newvalue .mesa 8 30 0.950000 0.000000 0.950000
xc_newvalue .mesa 8 31 0.950000 0.000000 0.950000
xc_newvalue .mesa 8 32 0.950000 0.000000 0.950000
xc_newvalue .mesa 8 33 0.950000 0.950000 0.000000
xc_newvalue .mesa 8 34 0.950000 0.950000 0.000000
xc_newvalue .mesa 8 35 0.950000 0.000000 0.000000
xc_newvalue .mesa 8 36 0.750000 0.750000 0.750000
xc_newvalue .mesa 8 37 0.000000 0.950000 0.950000
xc_newvalue .mesa 8 38 0.000000 0.950000 0.950000
xc_newvalue .mesa 8 39 0.750000 0.750000 0.750000
xc_newvalue .mesa 8 40 0.750000 0.750000 0.750000
xc_newvalue .mesa 8 41 0.750000 0.750000 0.750000
xc_newvalue .mesa 8 42 0.750000 0.750000 0.750000
xc_newvalue .mesa 8 43 0.750000 0.750000 0.750000
xc_newvalue .mesa 8 44 0.750000 0.750000 0.750000
xc_newvalue .mesa 8 45 0.750000 0.750000 0.750000
xc_newvalue .mesa 8 46 1.000000 1.000000 1.000000
xc_newvalue .mesa 8 47 1.000000 1.000000 1.000000
xc_newvalue .mesa 8 48 0.750000 0.750000 0.750000
xc_newvalue .mesa 8 49 0.750000 0.750000 0.750000
xc_newvalue .mesa 8 50 0.750000 0.750000 0.750000
xc_newvalue .mesa 8 51 0.750000 0.750000 0.750000
xc_newvalue .mesa 8 52 0.750000 0.750000 0.750000
xc_newvalue .mesa 8 53 0.750000 0.750000 0.750000
xc_newvalue .mesa 8 54 0.750000 0.750000 0.750000
xc_newvalue .mesa 8 55 0.000000 0.950000 0.950000
xc_newvalue .mesa 8 56 0.000000 0.950000 0.950000
xc_newvalue .mesa 8 57 0.750000 0.750000 0.750000
xc_newvalue .mesa 8 58 0.750000 0.750000 0.750000
xc_newvalue .mesa 8 59 0.750000 0.750000 0.750000
xc_newvalue .mesa 8 60 0.750000 0.750000 0.750000
xc_newvalue .mesa 8 61 0.750000 0.750000 0.750000
xc_newvalue .mesa 8 62 0.750000 0.750000 0.750000
xc_newvalue .mesa 8 63 0.750000 0.750000 0.750000
xc_newvalue .mesa 8 64 0.750000 0.750000 0.750000
xc_newvalue .mesa 8 65 0.750000 0.750000 0.750000
xc_newvalue .mesa 8 66 0.750000 0.750000 0.750000
xc_newvalue .mesa 8 67 0.750000 0.750000 0.750000
xc_newvalue .mesa 8 68 0.750000 0.750000 0.750000
xc_newvalue .mesa 8 69 0.750000 0.750000 0.750000
xc_newvalue .mesa 8 70 0.750000 0.750000 0.750000
xc_newvalue .mesa 8 71 0.750000 0.750000 0.750000
xc_newvalue .mesa 8 72 0.750000 0.750000 0.750000
xc_newvalue .mesa 8 73 0.750000 0.750000 0.750000
xc_newvalue .mesa 8 74 0.750000 0.750000 0.750000
xc_newvalue .mesa 8 75 0.750000 0.750000 0.750000
xc_newvalue .mesa 8 76 0.750000 0.750000 0.750000
xc_newvalue .mesa 8 77 0.750000 0.750000 0.750000
xc_newvalue .mesa 8 78 0.750000 0.750000 0.750000
xc_newvalue .mesa 8 79 1.000000 0.850000 0.000000
xc_newvalue .mesa 8 80 0.750000 0.750000 0.750000
xc_newvalue .mesa 8 81 0.750000 0.750000 0.750000
xc_newvalue .mesa 8 82 0.750000 0.750000 0.750000
xc_newvalue .mesa 8 83 0.750000 0.750000 0.750000
xc_newvalue .mesa 8 84 0.750000 0.750000 0.750000
xc_newvalue .mesa 8 85 0.750000 0.750000 0.750000
xc_newvalue .mesa 8 86 0.750000 0.750000 0.750000
xc_newvalue .mesa 8 87 0.750000 0.750000 0.750000
xc_newvalue .mesa 8 88 0.750000 0.750000 0.750000
xc_newvalue .mesa 8 89 0.750000 0.750000 0.750000
xc_newvalue .mesa 8 90 0.750000 0.750000 0.750000
xc_newvalue .mesa 8 91 0.750000 0.750000 0.750000
xc_newvalue .mesa 8 92 0.750000 0.750000 0.750000
xc_newvalue .mesa 8 93 0.750000 0.750000 0.750000
xc_newvalue .mesa 8 94 0.750000 0.750000 0.750000
xc_newvalue .mesa 8 95 0.750000 0.750000 0.750000
xc_newvalue .mesa 8 96 0.750000 0.750000 0.750000
xc_newvalue .mesa 8 97 0.750000 0.750000 0.750000
xc_newvalue .mesa 8 98 0.750000 0.750000 0.750000
xc_newvalue .mesa 8 99 0.931400 0.875500 0.801000
xc_newvalue .mesa 8 100 0.931400 0.875500 0.801000
xc_newvalue .mesa 26 1.000000 1.000000 1.000000 1.0
xc_newvalue .mesa 17 0.880000 1.000000 0.670000
# ------------------------------------------------------------------------
# END: Various colors
# ------------------------------------------------------------------------
xc_newvalue .mesa 0
# ------------------------------------------------------------------------
# BEGIN: Atomic radii
# ------------------------------------------------------------------------
xc_newvalue .mesa 4 0 0.532000
xc_newvalue .mesa 4 1 0.532000
xc_newvalue .mesa 4 2 0.532000
xc_newvalue .mesa 4 3 1.722000
xc_newvalue .mesa 4 4 1.246000
xc_newvalue .mesa 4 5 1.274000
xc_newvalue .mesa 4 6 1.078000
xc_newvalue .mesa 4 7 1.050000
xc_newvalue .mesa 4 8 1.022000
xc_newvalue .mesa 4 9 0.994000
xc_newvalue .mesa 4 10 0.994000
xc_newvalue .mesa 4 11 2.240000
xc_newvalue .mesa 4 12 1.960000
xc_newvalue .mesa 4 13 1.750000
xc_newvalue .mesa 4 14 1.554000
xc_newvalue .mesa 4 15 1.400000
xc_newvalue .mesa 4 16 1.456000
xc_newvalue .mesa 4 17 1.386000
xc_newvalue .mesa 4 18 1.372000
xc_newvalue .mesa 4 19 2.982000
xc_newvalue .mesa 4 20 2.436000
xc_newvalue .mesa 4 21 2.240000
xc_newvalue .mesa 4 22 1.960000
xc_newvalue .mesa 4 23 1.890000
xc_newvalue .mesa 4 24 1.960000
xc_newvalue .mesa 4 25 1.960000
xc_newvalue .mesa 4 26 1.960000
xc_newvalue .mesa 4 27 1.890000
xc_newvalue .mesa 4 28 1.890000
xc_newvalue .mesa 4 29 1.890000
xc_newvalue .mesa 4 30 1.890000
xc_newvalue .mesa 4 31 1.820000
xc_newvalue .mesa 4 32 1.750000
xc_newvalue .mesa 4 33 1.610000
xc_newvalue .mesa 4 34 1.610000
xc_newvalue .mesa 4 35 1.596000
xc_newvalue .mesa 4 36 1.568000
xc_newvalue .mesa 4 37 3.080000
xc_newvalue .mesa 4 38 2.800000
xc_newvalue .mesa 4 39 2.590000
xc_newvalue .mesa 4 40 2.170000
xc_newvalue .mesa 4 41 2.030000
xc_newvalue .mesa 4 42 2.030000
xc_newvalue .mesa 4 43 1.890000
xc_newvalue .mesa 4 44 1.820000
xc_newvalue .mesa 4 45 1.890000
xc_newvalue .mesa 4 46 1.960000
xc_newvalue .mesa 4 47 2.240000
xc_newvalue .mesa 4 48 2.170000
xc_newvalue .mesa 4 49 2.170000
xc_newvalue .mesa 4 50 1.974000
xc_newvalue .mesa 4 51 2.030000
xc_newvalue .mesa 4 52 1.960000
xc_newvalue .mesa 4 53 1.960000
xc_newvalue .mesa 4 54 1.834000
xc_newvalue .mesa 4 55 3.640000
xc_newvalue .mesa 4 56 2.800000
xc_newvalue .mesa 4 57 2.450000
xc_newvalue .mesa 4 58 2.170000
xc_newvalue .mesa 4 59 2.170000
xc_newvalue .mesa 4 60 2.170000
xc_newvalue .mesa 4 61 2.170000
xc_newvalue .mesa 4 62 2.170000
xc_newvalue .mesa 4 63 2.170000
xc_newvalue .mesa 4 64 2.170000
xc_newvalue .mesa 4 65 2.170000
xc_newvalue .mesa 4 66 2.170000
xc_newvalue .mesa 4 67 2.170000
xc_newvalue .mesa 4 68 2.170000
xc_newvalue .mesa 4 69 2.170000
xc_newvalue .mesa 4 70 2.170000
xc_newvalue .mesa 4 71 2.170000
xc_newvalue .mesa 4 72 2.170000
xc_newvalue .mesa 4 73 2.030000
xc_newvalue .mesa 4 74 1.890000
xc_newvalue .mesa 4 75 1.890000
xc_newvalue .mesa 4 76 1.820000
xc_newvalue .mesa 4 77 1.890000
xc_newvalue .mesa 4 78 1.890000
xc_newvalue .mesa 4 79 1.890000
xc_newvalue .mesa 4 80 2.100000
xc_newvalue .mesa 4 81 2.660000
xc_newvalue .mesa 4 82 2.520000
xc_newvalue .mesa 4 83 2.240000
xc_newvalue .mesa 4 84 2.170000
xc_newvalue .mesa 4 85 2.170000
xc_newvalue .mesa 4 86 2.170000
xc_newvalue .mesa 4 87 3.920000
xc_newvalue .mesa 4 88 2.016000
xc_newvalue .mesa 4 89 2.730000
xc_newvalue .mesa 4 90 2.170000
xc_newvalue .mesa 4 91 2.170000
xc_newvalue .mesa 4 92 2.170000
xc_newvalue .mesa 4 93 2.170000
xc_newvalue .mesa 4 94 2.170000
xc_newvalue .mesa 4 95 2.170000
xc_newvalue .mesa 4 96 2.170000
xc_newvalue .mesa 4 97 2.170000
xc_newvalue .mesa 4 98 2.170000
xc_newvalue .mesa 4 99 2.170000
xc_newvalue .mesa 4 100 2.170000
xc_newvalue .mesa 10004 0 0.000000
xc_newvalue .mesa 10004 1 0.399000
xc_newvalue .mesa 10004 2 0.399000
xc_newvalue .mesa 10004 3 1.291500
xc_newvalue .mesa 10004 4 0.934500
xc_newvalue .mesa 10004 5 0.955500
xc_newvalue .mesa 10004 6 0.808500
xc_newvalue .mesa 10004 7 0.787500
xc_newvalue .mesa 10004 8 0.766500
xc_newvalue .mesa 10004 9 0.745500
xc_newvalue .mesa 10004 10 0.745500
xc_newvalue .mesa 10004 11 1.680000
xc_newvalue .mesa 10004 12 1.470000
xc_newvalue .mesa 10004 13 1.312500
xc_newvalue .mesa 10004 14 1.165500
xc_newvalue .mesa 10004 15 1.050000
xc_newvalue .mesa 10004 16 1.092000
xc_newvalue .mesa 10004 17 1.039500
xc_newvalue .mesa 10004 18 1.029000
xc_newvalue .mesa 10004 19 2.236500
xc_newvalue .mesa 10004 20 1.827000
xc_newvalue .mesa 10004 21 1.680000
xc_newvalue .mesa 10004 22 1.470000
xc_newvalue .mesa 10004 23 1.417500
xc_newvalue .mesa 10004 24 1.470000
xc_newvalue .mesa 10004 25 1.470000
xc_newvalue .mesa 10004 26 1.470000
xc_newvalue .mesa 10004 27 1.417500
xc_newvalue .mesa 10004 28 1.417500
xc_newvalue .mesa 10004 29 1.417500
xc_newvalue .mesa 10004 30 1.417500
xc_newvalue .mesa 10004 31 1.365000
xc_newvalue .mesa 10004 32 1.312500
xc_newvalue .mesa 10004 33 1.207500
xc_newvalue .mesa 10004 34 1.207500
xc_newvalue .mesa 10004 35 1.197000
xc_newvalue .mesa 10004 36 1.176000
xc_newvalue .mesa 10004 37 2.310000
xc_newvalue .mesa 10004 38 2.100000
xc_newvalue .mesa 10004 39 1.942500
xc_newvalue .mesa 10004 40 1.627500
xc_newvalue .mesa 10004 41 1.522500
xc_newvalue .mesa 10004 42 1.522500
xc_newvalue .mesa 10004 43 1.417500
xc_newvalue .mesa 10004 44 1.365000
xc_newvalue .mesa 10004 45 1.417500
xc_newvalue .mesa 10004 46 1.470000
xc_newvalue .mesa 10004 47 1.680000
xc_newvalue .mesa 10004 48 1.627500
xc_newvalue .mesa 10004 49 1.627500
xc_newvalue .mesa 10004 50 1.480500
xc_newvalue .mesa 10004 51 1.522500
xc_newvalue .mesa 10004 52 1.470000
xc_newvalue .mesa 10004 53 1.470000
xc_newvalue .mesa 10004 54 1.375500
xc_newvalue .mesa 10004 55 2.730000
xc_newvalue .mesa 10004 56 2.100000
xc_newvalue .mesa 10004 57 1.837500
xc_newvalue .mesa 10004 58 1.627500
xc_newvalue .mesa 10004 59 1.627500
xc_newvalue .mesa 10004 60 1.627500
xc_newvalue .mesa 10004 61 1.627500
xc_newvalue .mesa 10004 62 1.627500
xc_newvalue .mesa 10004 63 1.627500
xc_newvalue .mesa 10004 64 1.627500
xc_newvalue .mesa 10004 65 1.627500
xc_newvalue .mesa 10004 66 1.627500
xc_newvalue .mesa 10004 67 1.627500
xc_newvalue .mesa 10004 68 1.627500
xc_newvalue .mesa 10004 69 1.627500
xc_newvalue .mesa 10004 70 1.627500
xc_newvalue .mesa 10004 71 1.627500
xc_newvalue .mesa 10004 72 1.627500
xc_newvalue .mesa 10004 73 1.522500
xc_newvalue .mesa 10004 74 1.417500
xc_newvalue .mesa 10004 75 1.417500
xc_newvalue .mesa 10004 76 1.365000
xc_newvalue .mesa 10004 77 1.417500
xc_newvalue .mesa 10004 78 1.417500
xc_newvalue .mesa 10004 79 1.417500
xc_newvalue .mesa 10004 80 1.575000
xc_newvalue .mesa 10004 81 1.995000
xc_newvalue .mesa 10004 82 1.890000
xc_newvalue .mesa 10004 83 1.680000
xc_newvalue .mesa 10004 84 1.627500
xc_newvalue .mesa 10004 85 1.627500
xc_newvalue .mesa 10004 86 1.627500
xc_newvalue .mesa 10004 87 2.940000
xc_newvalue .mesa 10004 88 1.512000
xc_newvalue .mesa 10004 89 2.047500
xc_newvalue .mesa 10004 90 1.627500
xc_newvalue .mesa 10004 91 1.627500
xc_newvalue .mesa 10004 92 1.627500
xc_newvalue .mesa 10004 93 1.627500
xc_newvalue .mesa 10004 94 1.627500
xc_newvalue .mesa 10004 95 1.627500
xc_newvalue .mesa 10004 96 1.627500
xc_newvalue .mesa 10004 97 1.627500
xc_newvalue .mesa 10004 98 1.627500
xc_newvalue .mesa 10004 99 1.627500
xc_newvalue .mesa 10004 100 1.627500
# ------------------------------------------------------------------------
# END: Atomic radii
# ------------------------------------------------------------------------
# ------------------------------------------------------------------------
# Various parameters
# ------------------------------------------------------------------------
xc_newvalue .mesa 13 1.050000000000000
xc_newvalue .mesa 6 0.400000000000000
xc_newvalue .mesa 7 0.600000000000000
xc_newvalue .mesa 9 1.000000000000000
xc_newvalue .mesa 10 1.000000000000000
xc_newvalue .mesa 18 1.000000000000000
xc_newvalue .mesa 10010 1.000000000000000
xc_newvalue .mesa 10009 1.000000000000000
xc_newvalue .mesa 11 6.000000000000000
xc_newvalue .mesa 19 0.100000000000000
xc_newvalue .mesa 24 15.000000000000000
xc_newvalue .mesa 28 3.000000000000000
xc_newvalue .mesa 29 2.500000000000000
xc_newvalue .mesa 10029 1.000000000000000
xc_newvalue .mesa 23 1.000000 1.000000 1.000000 1.000000
# ------------------------------------------------------------------------
# Various displays (i.e. checkbuttons of DISPLAY menu)
# ------------------------------------------------------------------------
array set check {pseudoDens 0 perspective 0 labels 0 depthcuing 0 crds 0 wigner 0 antialias 0 Hbonds 0 forces 0 frames 0 unibond 0 perpective 0}
CrdSist
AtomLabels
CrysFrames
Unibond
forceVectors .mesa
WignerSeitz
Perspective
# ------------------------------------------------------------------------
# Various displays (i.e. radiobuttons DISPLAY menu)
# ------------------------------------------------------------------------
array set radio {space {SpaceFill based on covalent radii} .mesa,bg #ffffff cellmode conv frames rods unitrep cell hexamode parapipedal ball {Balls based on covalent radii}}
xc_newvalue .mesa 2
# ------------------------------------------------------------------------
# load appropriate color-scheme
# ------------------------------------------------------------------------
array set colSh {slab_dir -z slabrange_min 0.00 slabrange_max 1.00 scheme atomic slab_colbas monocrome dist_r 1.0 dist_coltyp combined dist_alpha 0.65 dist_x 0.0 slab_fractional 1 slab_coltyp combined dist_y 0.0 dist_z 0.0 slab_alpha 0.65 dist_colbas monocrome}
ColorSchemeUpdate .mesa
# ------------------------------------------------------------------------
# BEGIN: Isosurface colors/transparency ...
# ------------------------------------------------------------------------
xc_setGLparam material -what isosurf pos front \
-shininess 0.000000 \
-specular {0.000000 0.000000 0.000000 0.000000} \
-ambient {0.000000 0.000000 0.000000 0.000000} \
-diffuse {0.000000 0.000000 0.000000 0.000000} \
-emission {0.000000 0.000000 0.000000 0.000000}
xc_setGLparam material -what isosurf pos back \
-shininess 0.000000 \
-specular {0.000000 0.000000 0.000000 0.000000} \
-ambient {0.000000 0.000000 0.000000 0.000000} \
-diffuse {0.000000 0.000000 0.000000 0.000000} \
-emission {0.000000 0.000000 0.000000 0.000000}
xc_setGLparam material -what isosurf neg front \
-shininess 0.000000 \
-specular {0.000000 0.000000 0.000000 0.000000} \
-ambient {0.000000 0.000000 0.000000 0.000000} \
-diffuse {0.000000 0.000000 0.000000 0.000000} \
-emission {0.000000 0.000000 0.000000 0.000000}
xc_setGLparam material -what isosurf neg back \
-shininess 0.000000 \
-specular {0.000000 0.000000 0.000000 0.000000} \
-ambient {0.000000 0.000000 0.000000 0.000000} \
-diffuse {0.000000 0.000000 0.000000 0.000000} \
-emission {0.000000 0.000000 0.000000 0.000000}
xc_setGLparam material -what isosurf one front \
-shininess 0.000000 \
-specular {0.000000 0.000000 0.000000 0.000000} \
-ambient {0.000000 0.000000 0.000000 0.000000} \
-diffuse {0.000000 0.000000 0.000000 0.000000} \
-emission {0.000000 0.000000 0.000000 0.000000}
xc_setGLparam material -what isosurf one back \
-shininess 0.000000 \
-specular {0.000000 0.000000 0.000000 0.000000} \
-ambient {0.000000 0.000000 0.000000 0.000000} \
-diffuse {0.000000 0.000000 0.000000 0.000000} \
-emission {0.000000 0.000000 0.000000 0.000000}
xc_setGLparam blendfunc -what isosurf -sfunc GL_SRC_ALPHA -dfunc GL_ONE_MINUS_SRC_ALPHA
# ------------------------------------------------------------------------
# END: Isosurface colors/transparency ...
# ------------------------------------------------------------------------
# ------------------------------------------------------------------------
# rotation matrix and zooming factor, and translation displacements
# ------------------------------------------------------------------------
xc_rotationmatrix set 1.000000000000000e+00 0.000000000000000e+00 0.000000000000000e+00 0.000000000000000e+00 0.000000000000000e+00 1.000000000000000e+00 0.000000000000000e+00 0.000000000000000e+00 0.000000000000000e+00 0.000000000000000e+00 1.000000000000000e+00 0.000000000000000e+00 0.000000000000000e+00 0.000000000000000e+00 0.000000000000000e+00 1.000000000000000e+00
xc_translparam set 0.000000000000000e+00 0.000000000000000e+00 2.038778501368373e+00
# this is used to force the update of display
.mesa cry_toglzoom 0.0
# ------------------------------------------------------------------------
# Anti-Aliasing & Depth-Cuing & PseudoDensity (these are time consuming)
# ------------------------------------------------------------------------
DepthCuing; PseudoDensity; AntiAlias
# ------------------------------------------------------------------------
# reset cursor
# ------------------------------------------------------------------------
set xcCursor(dont_update) 0
ResetCursor
destroy \$wait_window
# now print the MO
update
scripting::printToFile $i-$j.png
exit
EOF
# end of embedded xcrysden.script
# now we load the just created xcrysden.script into xcrysden
xcrysden -s xcrysden.script
done
done
rm xcrysden.script
mkdir ini-fin
mv conv-ini-fin.jpg ini-fin
mv *.png ini-fin
mv result.lipai ini-fin
mogrify -trim ini-fin/*.png
mogrify -trim ini-fin/*.jpg
plotdist.sh¶
This bash script supervises distant changes of structures in VASP optimization.
#!/bin/bash
touch p1 p2;
touch dist.conv
Num=`awk 'NR==7{for(i=1;i<=NF;i++) a=$i+a;}END{print a}' XDATCAR`
Lnum=`wc XDATCAR|awk '{print $1}'`
((n=(Lnum-7)/(Num+1)))
head -7 XDATCAR >p1
awk -v num="$Num" 'NR==9,NR==(num+1)+7{print $0}' XDATCAR >>p1
for((i=1;i<n;i++))
do
cat p1 >p2
head -7 XDATCAR >p1
((n1=9+(Num+1)*i))
((n2=(i+1)*(Num+1)+7))
sed -n ''$n1','$n2'p' XDATCAR >>p1
echo -e $i"\t"`dist.pl p1 p2 ` >>dist.conv
done
# plot
gnuplot <<EOF
set grid
set term post
set output 'b.ps'
set xlabel 'Ion steps'
set title 'Distance between each ion steps'
unset key
set ylabel 'dist(Angst)'
plot 'dist.conv' u 1:2 w l lw 2 lc rgb "blue"
EOF
gs -sDEVICE=jpeg -r300 -sPAPERSIZE=a4 -dBATCH -dNOPAUSE -sOutputFile=dist.jpg b.ps
convert -rotate 90 dist.jpg dist.jpg
mogrify -trim dist.jpg
gnuplot <<EOF
set term dumb
plot 'dist.conv' w l t "Dist "
EOF
rm b.ps p1 p2 dist.conv
plot Band Structure¶
The bash script extracts data of band structure from OUTCAR file, and write the output file band.dat The python script use the band.dat to plot band structure with the output file band.jpg
band.sh
#edit by lipai@mail.ustc.edu.cn
#get num of band
num=`grep NBANDS OUTCAR |awk '{print $NF}'`
E_fermi=`grep E-fermi OUTCAR | awk '{print $3}'`
grep "band No." OUTCAR -A $num >temp
echo "" >>temp
split -l 26 -d temp bands
for i in bands*; do sed -i '$d' $i;awk -v num=$i 'BEGIN{printf("%s\t",num)};{if(NR>1) printf("%f\t",$2)};END{printf("\n")}' $i >>band.dat; done
sed -i s/bands// band.dat
rm bands* temp
python plot_band band.dat $num $E_fermi
plot_band.py
import sys
import matplotlib as mpl
import matplotlib.pyplot as plt
import numpy as np
mpl.rcParams['agg.path.chunksize']=10000
arg_num=len(sys.argv)
if arg_num!=5:
print("wrong inputs!!!\n")
print("format:")
print("plot_band filename band_num E_fermi")
exit()
filename=sys.argv[1]
band_num=int(sys.argv[2])
E_fermi=float(sys.argv[3])
split=int(sys.argv[4])
data=np.loadtxt(filename)
#plt.figure(figsize=(10,10))
for i in range(1,band_num+1):
#plt.plot(data[:,0],data[:,i],'.b')
plt.plot(data[:,0],data[:,i]-E_fermi,'-g',linewidth=3,alpha=0.4)
min=10
max=-10
for (x,y),val in np.ndenumerate(data):
val=val-E_fermi
if val>0 and val<min:
min=val
minx=x
miny=y
elif val<0 and val>max:
max=val
maxx=x
maxy=y
plt.plot(minx,min,'or')
plt.plot(maxx,max,'or')
ax=plt.gca()
ax.annotate(" "+str(minx)+" "+str(min),xy=(minx+x/30,min),color='r',bbox=dict(boxstyle='round,pad=0.2', fc='w', alpha=0.7))
ax.annotate(" "+str(maxx)+" "+str(max),xy=(maxx+x/30,max),color='r',bbox=dict(boxstyle='round,pad=0.2', fc='w', alpha=0.7))
if minx==maxx:
text="Direct band gap with gap of "+str(min-max)+" eV"
else:
text="Indirect band gap with gap of "+str(min-max)+" eV"
plt.title(text)
plt.ylabel('Energy/eV')
plt.ylim(-4,4)
plt.xlim(data[0,0],data[-1,0])
plt.axhline(0.0,label='E-fermi',color='r',alpha=0.2)
for i in range(1,split):
plt.axvline(i*data[-1,0]/split,color='k',alpha=0.2)
#plt.savefig("band.jpg",dpi=300)
plt.show()
calculate thermodymamics¶
This python script calculate ZPE/TS/H using frequency info in OUTCAR
from monty.re import regrep
def get_frequency(outcar):
pattern = {"frequency": "^\s+[\d]+\s+f+\s+=\s+([\d\-\.]+)\sTHz+\s+([\d\-\.]+)\s2PiTHz+\s+([\d\-\.]+)\s+cm-1+\s+([\d\-\.]+)+\smeV"}
d = regrep(outcar, pattern)
# float string and rename index
data = []
for index, d in enumerate(d['frequency']):
d[0] = [float(i) for i in d[0]]
d[1] = '{} f'.format(index + 1)
data.append(d)
return data
dd = get_frequency('OUTCAR')
meV_data = [d[-1] for d, index in dd]
l=len(meV_data)
Kb = 8.6173324E-2 #Boltzmann constant unit is meV/K
T = 298 #Absolute temperature unit is K
R = 8.3144598 #Gas constant unit is J/(k.mol)
Tran = 1.0364E-5 #J/mol unit to eV
D = []
S = []
C = []
import numpy as np
for i in xrange(0,l):
d =meV_data[i]/(Kb*T)
D.append(d)
s = (d/(np.exp(d)-1))-np.log1p(-np.exp(-d))
S.append(s)
c = meV_data[i]/(Kb*(np.exp(d)-1))
C.append(c)
entropy= sum(S)
TS= R*T*entropy*Tran
ZPE= 0.5*sum(meV_data)*1E-3
U = R*Tran*sum(C)
print "ZPE =", ZPE, "eV"
print "TS =", TS, "eV"
print "U =", U, "eV"
fix_atoms¶
#!/usr/bin/python
#The script is written by liren,Email: llr@mail.ustc.edu.cn
import sys
aname=sys.argv[1]
bname=sys.argv[2]
c=float(sys.argv[3])
a=open(aname)
b=open(bname,'w')
#line=a.readlines()
title=a.readline()
whatever=a.readline()
latt1=a.readline()
latt2=a.readline()
latt3=a.readline()
atom=a.readline()
natom=a.readline()
direct=a.readline()
b.write('%s\n'%title.strip())
b.write('%s\n'%whatever.strip())
b.write('%s\n'%latt1.strip())
b.write('%s\n'%latt2.strip())
b.write('%s\n'%latt3.strip())
b.write('%s\n'%atom.strip())
b.write('%s\n'%natom.strip())
b.write('%s\n'%'Selective dynamics')
b.write('%s\n'%direct.strip())
for line in a.readlines():
num=line.split()
if float(num[2])<=c:
b.write('%-16s%-16s%-16s%3s%3s%3s\n'%(num[0],num[1],num[2],'F','F','F'))
else:
b.write('%-16s%-16s%-16s%3s%3s%3s\n'%(num[0],num[1],num[2],'T','T','T'))
b.close()
Idpp method for creating initial NEB path¶
IDPP method is introduced by Hannes Jónsson in his paper “Improved initial guess for minimum energy path calculations”
#!/usr/bin/python
#lipai@mail.ustc.edu.cn
import numpy as np
import os
import time
step_init=0.0001
images=input("input num of images: ")
ininame=raw_input("ini structure: ")
finname=raw_input("fin structure: ")
#read ini structure
fileopen=open(ininame,'r')
ini_data=fileopen.readlines()
head=ini_data[:9]
atom_num=sum(map(int,head[6].split()))
ini_data=ini_data[9:9+atom_num]
tmp=[]
fix=[]
for i in range(atom_num):
tmp.append(map(float,ini_data[i].split()[0:3]))
if(ini_data[i].split()[4]=='F'): fix.append(1)
else: fix.append(0)
pos_a=np.array(tmp)
#read fin stucture
fileopen=open(finname,'r')
fin_data=fileopen.readlines()
fin_data=fin_data[9:9+atom_num]
tmp=[]
for i in fin_data:
tmp.append(map(float,i.split()[0:3]))
pos_b=np.array(tmp)
#correction of periodic boundary condition
for i in range(atom_num):
for j in range(3):
if(pos_a[i][j]-pos_b[i][j]>0.5): pos_a[i][j]-=1
if(pos_a[i][j]-pos_b[i][j]<-0.5): pos_b[i][j]-=1
#get dist matrix and linear interpolation
dist_a=np.zeros([atom_num,atom_num])
dist_b=np.zeros([atom_num,atom_num])
for i in range(atom_num):
for j in range(atom_num):
tmp_a=0
tmp_b=0
for k in range(3):
tmp_a+=(pos_a[i][k]-pos_a[j][k])**2
tmp_b+=(pos_b[i][k]-pos_b[j][k])**2
dist_a[i,j]=np.sqrt(tmp_a)
dist_b[i,j]=np.sqrt(tmp_b)
dist_im=np.zeros([images,atom_num,atom_num])
pos_im=np.zeros([images,atom_num,3])
for i in range(images):
dist_im[i]=dist_a+(i+1.0)*(dist_b-dist_a)/(images+1.0)
pos_im[i]=pos_a+(i+1.0)*(pos_b-pos_a)/(images+1.0)
#optimization using steepest descent method
pos_tmp=np.zeros([atom_num,3])
dist_tmp=np.zeros([atom_num,atom_num])
s0=np.zeros(images)
s1=np.zeros(images)
flag=np.zeros(images)
for im in range(images):
if(flag[im]==1): continue
step=step_init
print "generate image "+str(im+1)
loop=0
while(1):
for i in range(atom_num): #get dist_tmp
for j in range(atom_num):
if(i==j):
dist_tmp[i,j]=10
else:
tmp=0
for k in range(3):
tmp+=(pos_im[im][i][k]-pos_im[im][j][k])**2
dist_tmp[i,j]=np.sqrt(tmp)
for i in range(atom_num):
for sigma in range(3):
grad=0
for j in range(atom_num):
if(j!=i):
grad+=2*(dist_im[im][i][j]-dist_tmp[i][j])*(pos_im[im][i][sigma]-pos_im[im][j][sigma])\
*(2*dist_im[im][i][j]-dist_tmp[i][j])/dist_tmp[i,j]**5
pos_tmp[i][sigma]=pos_im[im][i][sigma]+step*grad
pos_im[im]=pos_tmp
#judge convergence
s0[im]=s1[im]
s1[im]=0
for i in range(atom_num):
for j in range(i):
s1[im]+=(dist_im[im][i][j]-dist_tmp[i][j])**2/dist_tmp[i][j]**4
loop+=1
print "loop:%d"%loop
if(abs(s0[im]-s1[im])<0.01):
print "Converged!"
flag[im]=1
break
if(loop>1 and s1[im]>s0[im]): step=step/3
#mkdir and generate poscar files
if(images+1<10): num='0'+str(images+1)
else: num=str(images+1)
os.system("mkdir 00")
os.system("cp p1 00/POSCAR")
os.system("mkdir "+num)
os.system("cp p2 "+num+"/POSCAR")
for i in range(images):
if(i+1<10): num='0'+str(i+1)
else: num=str(i+1)
os.system("mkdir "+num)
data=pos_im[i].tolist()
filename=num+"/POSCAR"
f=file(filename,"a+")
f.writelines(head)
for j in range(atom_num):
line=map(str,data[j])
line=" ".join(line)
if(fix[j]==1): line=line+' F F F\n'
else: line=line+' T T T\n'
f.write(line)
f.close()
plot DOS & Band Structure¶
This python script plots dos and band using VASP results
画DOS
%matplotlib inline
from pymatgen.io.vasp import Vasprun
from pymatgen.electronic_structure.plotter import DosPlotter
v = Vasprun('Si-dos/vasprun.xml')
tdos = v.tdos
plotter = DosPlotter()
plotter.add_dos("Total DOS", tdos)
plotter.show(xlim=[-5, 5], ylim=[0, 4])
画pDOS
%matplotlib inline
from pymatgen.io.vasp import Vasprun
from pymatgen.electronic_structure.plotter import DosPlotter
v = Vasprun('Si-dos/vasprun.xml')
cdos = v.complete_dos
element_dos = cdos.get_element_dos()
plotter = DosPlotter()
plotter.add_dos_dict(element_dos)
plotter.show(xlim=[-5, 5], ylim=[0, 1])
画band structure
from pymatgen.io.vasp import Vasprun, BSVasprun
from pymatgen.electronic_structure.plotter import BSPlotter
v = BSVasprun("Si-band/vasprun.xml")
bs = v.get_band_structure(kpoints_filename="Si-band/KPOINTS",line_mode=True)
plt = BSPlotter(bs)
plt.get_plot(vbm_cbm_marker=True,ylim=(-3,3))
plt.show()
Warning
gnuplot should be installed for this script, otherwise, it can not plot.
Scripts for job management system on super-computer¶
PBS system¶
grep lipai |cut -d . -f 1|xargs qstat -f |egrep "state|DIR"|cut -d K -f 2|cut -d , -f 1|cut -d = -f 2
LSF system¶
echo ================================================================================================
bjobs
echo ================================================================================================
echo
echo ================================================================================================
bjobs -l |egrep -A 1 "CWD"|sed -n 'N;N;s/\n//g'p|awk '{if(NR%2==1) print $0}'|sed 's/ //g'|cut -d , -f 2|cut -d "<" -f 2|cut -d ">" -f 1
echo ================================================================================================
echo
echo ================================================================================================
#bjobs -l |egrep -A 1 "CWD"|xargs -l3|cut -d , -f 2|awk '{if($1=="Execution") {print "Run";print ""} else print $0}'\
#|sed 's/ //g'|cut -d "<" -f 2|cut -d ">" -f 1
bqueues|grep zyli
echo ================================================================================================
bhosts |grep -A 8 node291
bhosts |grep node300
echo ================================================================================================
for i in `seq 291 300`
do
lsload |grep node$i
done
echo ================================================================================================
YH system¶
yhqueue |grep yang
yhcontrol show jobs|grep -B 13 "yangjl"|grep Command|cut -d = -f 2
kMC simulation¶
Scripts¶
Framework of my kMC¶
This is a framework of lipai’s kMC code
#include<stdio.h>
#include<math.h>
#include<stdlib.h>
#include<vector>
#include<time.h>
#define N_limit 80000000
using namespace std;
const int multiplier=50;
const int step_out=10000000;
const int row=2000;
const int col=2000;
const int c_num=10;
double t=0.0;
//const int G_row=200;
//const int G_col=200;
const double kT=0.112;
int num_evt[2]={0};
int mesh[row][col][2]={0};// save species info
int species_num[2]={0};
int num_diff[2]={0}; //save diffu events
double prob_each[3][4]={
//0 1 2
//diff C C2 e1
/* {0.00,0.00,0.00,0.00},
{0.20,1.00,0.00,0.00},//C
{0.00,0.00,0.00,2.30},//C2 e1:C2-> C+C
*/
{0.00, 0.00, 0.00, 0.00},
{0.50, 0.90, 0.00, 0.00},
{0.49, 0.00, 0.00, 2.80},
};
// C2->C+C 2.7eV, C+C->C2 ???
int sp_merge[2][2]={\
//1 2
//C C2
{2,0},//C
{0,0},//C2
};
int sp_e[2][2]={\
{0,0}, //C
{1,1}}; //C2
int evt_emerge[2][2]={\
//1 2
//C C2
{1,-1},//C
{-1,-1},//C2
};
int evt_spec[2]={-1,2};
int p[3][2];
int q[3][2];
void err(char* erro);
void init();
int get_nb(int x,int y,int xy,int i);
void find_evt();
void diffusion(int sp,int s_n,int evt);
void merge(int sp1,int s_n,int evt);
void spec_del(int sp,int s_n);
void creat(int sp,int x,int y);
void spec_e(int sp,int s_n,int evt);
void update_evt_nb(int x,int y);
void statistic(int n,int N);
int get_sp(int s_x,int s_y);
void add_c(int sp);
class species{
public:
species(int sp,int x,int y){
spec=sp;
s_x=x;
s_y=y;
prob[3]=prob_each[spec][3]; //desorp or decomp events
check_event();
}
void check_event(){
double r9;
for(int i=0;i<3;i++){ //3 directions
nb[i][0]=get_nb(s_x,s_y,0,i);
nb[i][1]=get_nb(s_x,s_y,1,i);
event[i]=get_sp(nb[i][0],nb[i][1]);
// if(event[i]==100) prob[i]=0; //2017 added for h
// else prob[i]=prob_each[spec][event[i]];
prob[i]=prob_each[spec][event[i]];
}
prob_s=0;
for(int i=0;i<4;i++) prob_s+=prob[i];
}
int update(){
int i;
double evt_sum=0;
double evt_add=0;
double evt_rand;
for(i=0;i<4;i++)
if(prob[i]!=0) evt_sum+=prob[i];
evt_rand=evt_sum*(double)rand()/RAND_MAX;
for(i=0;i<4;i++){
if(prob[i]!=0) evt_add+=prob[i];
if(evt_add>=evt_rand){
return i;
}
}
err("spec.update()");
return 0;
}
double prob[4],prob_s;
int event[3],nb[3][2];
int spec,s_x,s_y;
};
vector<vector<species> > spec_list(2);
int main(){
init();
int init_h=0,h;
for(int i=0;i<2;i++) init_h+=species_num[i];
for(int m=0;m<multiplier;m++){
for(int i=0;i<N_limit;i++){
// h=init_h;
find_evt();
// init_h=0;
// for(int j=0;j<2;j++) init_h+=species_num[j];
//if(init_h!=h) statistic(i,m);
if(i%step_out==0) statistic(i,0);
}
}
statistic(-1,0);
for(int i=0;i<2;i++){
printf("num_evt %d:\t%d\n",i,num_evt[i]);
printf("num_diff %d:\t%d\n",i,num_diff[i]);
}
return 0;
}
void err(char* erro){
printf("%s error!!!\n",erro);
scanf("%s",&erro);
}
void init(){
int i,j;
srand((unsigned)time(NULL));
for(i=1;i<3;i++) //2 species
for(j=0;j<4;j++)
if(prob_each[i][j]!=0)
prob_each[i][j]=2.71*pow(10,13)*exp(-prob_each[i][j]/kT);
for(i=1;i<3;i++)
for(j=1;j<3;j++)
if(prob_each[i][j]!=0)
prob_each[i][j]=prob_each[i][j]/2;
p[0][0]=-1;p[0][1]=1;// p for x even
p[1][0]=-1;p[1][1]=0;
p[2][0]=1;p[2][1]=0;
q[0][0]=-1;q[0][1]=0;//q for x odd
q[1][0]=1;q[1][1]=0;
q[2][0]=1;q[2][1]=-1;
for(i=0;i<c_num;i++) add_c(1);
//for(int i=row/2-G_row/2-1;i<row/2+G_row/2;i++)
//for(int j=col/2-G_col/2-1;j<col/2+G_col/2;j++)
//mesh[i][j][0]=9;
//mesh[row/2-G_row/2-1][col/2-G_col/2-1][0]=100;
//mesh[row/2-G_row/2-1][col/2+G_col/2][0]=100;
//mesh[row/2+G_row/2][col/2-G_col/2-1][0]=100;
//mesh[row/2+G_row/2][col/2+G_col/2][0]=100;
}
int get_nb(int x,int y,int xy,int i){
int x_temp,y_temp;
if(xy==0){
if(x%2==0) x_temp=x+p[i][0];
else x_temp=x+q[i][0];
if(x_temp>=row) return x_temp-=row;
else if(x_temp<0) return x_temp+=row;
else return x_temp;
}
else {
if(x%2==0) y_temp=y+p[i][1];
else y_temp=y+q[i][1];
if(y_temp>=col) return y_temp-=col;
else if(y_temp<0) return y_temp+=col;
else return y_temp;
}
}
void find_evt(){
double prob_rand;
double prob_sum=0;
double prob_add=0;
int drct,s_n;
double delta_t;
for(int i=0;i<2;i++)
if(!spec_list[i].empty())
for(int j=0;j<spec_list[i].size();j++) prob_sum+=spec_list[i][j].prob_s;
//printf("prob_s:\t%f\n",prob_sum);
//
do{
delta_t=(double)rand()/RAND_MAX;
}while(delta_t==0);
delta_t=-log(delta_t)/prob_sum;
t+=delta_t;
prob_rand=prob_sum*(double)rand()/RAND_MAX;
for(int i=0;i<2;i++)
if(!spec_list[i].empty()){
for(int j=0;j<spec_list[i].size();j++){
prob_add+=spec_list[i][j].prob_s;
if(prob_add>=prob_rand){
drct=spec_list[i][j].update();
if(drct<3){
if(spec_list[i][j].event[drct]==0){
diffusion(i+1,j,drct);
num_diff[i]++;
//printf("diffuse cxhy");
//
//printf("get it1\n");
}
else if(spec_list[i][j].event[drct]<=2){
merge(i+1,j,drct);
//printf("get it2\n");
//printf("merge cxhy\n");
}
else if(spec_list[i][j].event[drct]==100)
return;
else err("error in update1");
}
else {
spec_e(i+1,j,drct);
if(evt_spec[i]!=-1) num_evt[evt_spec[i]-1]++;
}
return;
}
}
}
//err("prob");
}
void diffusion(int sp,int s_n,int drct){
int x1,y1,x2,y2;
x1=spec_list[sp-1][s_n].s_x;
y1=spec_list[sp-1][s_n].s_y;
x2=spec_list[sp-1][s_n].nb[drct][0];
y2=spec_list[sp-1][s_n].nb[drct][1];
if(mesh[x1][y1][0]!=sp||mesh[x2][y2][0]!=0){ //for debug
//printf("in diffu,sp,x1,y1,x2,y2: %d\t%d\t%d\t%d\t%d\n",sp,x1,y1,x2,y2);
//if(mesh[x1][y1][0]!=sp) {printf("sp=%d,mesh=%d\n",sp,mesh[x1][y1][0]); err("diffusion1");}
if(mesh[x1][y1][0]!=sp) printf("dif1");
if(mesh[x2][y2][0]!=0) printf("dif1");
printf("drct=%d\tsp1=%d,mesh[0]=%d\n",drct,sp,mesh[x1][y1][0]);
printf("event[drct]:%d",spec_list[sp-1][s_n].event[drct]);
printf("sp2=%d,mesh2[1]=%d\n",mesh[x2][y2][0],mesh[x2][y2][1]);
printf("x2=%d,y2[1]=%d\n",spec_list[sp-1][mesh[x2][y2][1]].s_x,spec_list[sp-1][mesh[x2][y2][1]].s_y);
for(int i=0;i<3;i++){
printf("evt[i]=%d\n",spec_list[sp-1][s_n].event[i]);
update_evt_nb(x2,y2);
printf("evt[i]=%d\n",spec_list[sp-1][s_n].event[i]);
spec_list[sp-1][s_n].check_event();
printf("evt[i]=%d\n",spec_list[sp-1][s_n].event[i]);
printf("neigbour of 1: %d\t%d\n",get_nb(x1,y1,0,i),get_nb(x1,y1,1,i));
printf("neigbour of 2: %d\t%d\n",get_nb(x2,y2,0,i),get_nb(x2,y2,1,i));
}
err("diffusion");
//}
err("diffusion");
}
mesh[x2][y2][0]=mesh[x1][y1][0];
mesh[x2][y2][1]=mesh[x1][y1][1];
mesh[x1][y1][0]=0;
mesh[x1][y1][1]=0;
spec_list[sp-1][s_n].s_x=x2;
spec_list[sp-1][s_n].s_y=y2;
spec_list[sp-1][s_n].check_event();
update_evt_nb(x1,y1);
update_evt_nb(x2,y2);
}
void merge(int sp1,int s_n,int drct){
int sp2,sp_crt;
int x1,y1,x2,y2;
x1=spec_list[sp1-1][s_n].s_x;
y1=spec_list[sp1-1][s_n].s_y;
x2=spec_list[sp1-1][s_n].nb[drct][0];
y2=spec_list[sp1-1][s_n].nb[drct][1];
sp2=spec_list[sp1-1][s_n].event[drct];
spec_del(sp1,s_n);
spec_del(sp2,mesh[x2][y2][1]);
sp_crt=sp_merge[sp1-1][sp2-1];
if(evt_emerge[sp1-1][sp2-1]!=-1) num_evt[evt_emerge[sp1-1][sp2-1]-1]++;
if(sp_crt>0) creat(sp_crt,x1,y1);
else if(sp_crt==0) err("merge");
}
void spec_del(int sp,int s_n){
int x,y;
int spec_end;
species_num[sp-1]--;
x=spec_list[sp-1][s_n].s_x;
y=spec_list[sp-1][s_n].s_y;
if(sp!=mesh[x][y][0]) err("spec_del1");
if(s_n!=mesh[x][y][1]) err("spec_del2");
spec_end=spec_list[sp-1].size()-1;
if(s_n!=spec_end){
//mesh[spec_list[sp-1][spec_end].s_x][spec_list[sp-1][spec_end].s_y][0]=sp;
mesh[spec_list[sp-1][spec_end].s_x][spec_list[sp-1][spec_end].s_y][1]=s_n;
swap(spec_list[sp-1][s_n],spec_list[sp-1][spec_end]);
}
spec_list[sp-1].pop_back();
mesh[x][y][0]=0;
mesh[x][y][1]=0;
update_evt_nb(x,y);
}
void creat(int sp,int x,int y){
species_num[sp-1]++;
mesh[x][y][0]=sp;
mesh[x][y][1]=spec_list[sp-1].size();
spec_list[sp-1].push_back(species(sp,x,y));
update_evt_nb(x,y);
}
void spec_e(int sp,int s_n,int drct){
int sp1,sp2,xn,yn;
int x,y;
x=spec_list[sp-1][s_n].s_x;
y=spec_list[sp-1][s_n].s_y;
sp1=sp_e[sp-1][0];
sp2=sp_e[sp-1][1];
// revised 0925
if(sp1==0) err("spec_e");
spec_del(sp,s_n);
creat(sp1,x,y);
for(int i=0;i<3;i++){
xn=get_nb(x,y,0,i);
yn=get_nb(x,y,1,i);
if(mesh[xn][yn][0]==0){
creat(sp2,xn,yn);
break;
}
}
}
void update_evt_nb(int x,int y){
int sp,s_n;
int x_temp,y_temp;
for(int i=0;i<3;i++){
x_temp=get_nb(x,y,0,i);
y_temp=get_nb(x,y,1,i);
sp=mesh[x_temp][y_temp][0];
if(sp>0&&sp<=8){
s_n=mesh[x_temp][y_temp][1];
spec_list[sp-1][s_n].check_event();
}
}
}
void statistic(int n,int N){
int nn[2];
for(int i=0;i<2;i++) nn[i]=spec_list[i].size();
printf("time/c/c2:\t%e\t%d\t%d\n",t,species_num[0],species_num[1]);
}
int get_sp(int x,int y){
return mesh[x][y][0];
}
void add_c(int sp){
int x,y,xn,yn;
int flag=0;
do{
do{
x=(int)(row*(double)rand()/RAND_MAX);
y=(int)(col*(double)rand()/RAND_MAX);
}while(mesh[x][y][0]!=0);
for(int i=0;i<3;i++){
xn=get_nb(x,y,0,i);
yn=get_nb(x,y,1,i);
if(mesh[xn][yn][0]==0){
flag=1;
break;
}
}
}while(flag==0);
creat(sp,x,y);
}
Machine Learning tools¶
Scripts¶
Extract data from Gaussian output file *.log to series XSF files¶
grep "SCF Done" *.log|awk '{print $5}' >energy.dat
cat *.log |sed -n '/concise/,/Rotational/p'| awk '{if($1 ~ /^[0-9]+$/) print $0;}' >stru.dat
natoms=`grep NAtoms *.log |head -1 |awk '{print $2}'`
split -l $natoms -a 3 -d stru.dat g
mkdir str
rm stru.dat
mv g* energy.dat str
cd str
mv g09.pbs ..
sed -i '$d' energy.dat
ls g* >files
a=`tail -1 files`
rm $a
a=`tail -2 files|head -1 `
rm $a
rm files
a=1
for i in g*
do
e=`head -n $a energy.dat |tail -1`
echo "# total energy = $e a.u." >$i.xsf
echo " " >>$i.xsf
echo "ATOMS" >>$i.xsf
sort -nk 2 $i |awk '{if($2=="1") printf("H %f %f %f 0 0 0\n",$4,$5,$6);
else printf("O %f %f %f 0 0 0\n",$4,$5,$6)} ' >>$i.xsf
a=$(($a+1))
done
cd ..
mkdir xsf
mv str/*.xsf xsf
Convert lammps structure to XYZ file¶
#!/bin/sh
#lipai@mail.ustc.edu.cn
if [ -d "xyz" ] ; then rm -rf xyz; fi
if [ -f "movie.xyz" ]; then rm movie.xyz;fi
for i in `seq 0 4000`
do
if [ -f "POSCAR-$i" ];then
pos2xyz.pl POSCAR-$i
cat POSCAR-$i.xyz >>movie.xyz
else
break
fi
done
mkdir xyz
mv *.xyz xyz
Extract data from OUTCAR to series XSF files¶
#!/bin/sh
#lipai@mail.ustc.edu.cn
#creat xsf files using OUTCAR and POSCAR
if [ $# = 0 ]; then
out="OUTCAR"
else
out=$1
fi
echo $out
rm *.temp
awk '/POSITION/,/drift/{
if(NF==6) print $0
}' $out > pos.temp
grep "energy without " $out |awk '{print $4}' >energy.temp
#head -5 POSCAR |tail -3 > primvec.temp
#atom_1=`awk '{if(NR==6) print $1}' POSCAR`
#atom_2=`awk '{if(NR==6) print $2}' POSCAR`
#num_1=`awk '{if(NR==7) print $1}' POSCAR`
#num_2=`awk '{if(NR==7) print $2}' POSCAR`
grep -A 3 "direct lattice vectors" $out \
|head -4|tail -3 |awk '{printf("%f %f %f\n",$1,$2,$3)}' >primvec.temp
atom_1=`grep "POTCAR" $out |awk 'NR==1{print $3}'`
atom_2=`grep "POTCAR" $out |awk 'NR==2{print $3}'`
num_1=`grep "ions per type" $out |head -1 \
|awk '{print $5}' `
num_2=`grep "ions per type" $out |head -1 \
|awk '{print $6}' `
num_atom=`echo "$num_1+$num_2" |bc`
lines=`wc pos.temp|awk '{print $1}'`
num_str=`echo "$lines/$num_atom" |bc`
for i in `seq $num_1`
do
echo "$atom_1" >> type.temp
done
for i in `seq $num_2`
do
echo "$atom_2" >> type.temp
done
echo $atom_1 $num_1
echo $atom_2 $num_2
echo "all $num_atom"
echo "num of str: $num_str"
for i in `seq $num_str`
do
energy=`head -n $i energy.temp|tail -1`
echo "# total energy = $energy eV" >> str_$i.xsf
echo " " >> str_$i.xsf
echo "CRYSTAL" >> str_$i.xsf
echo "PRIMVEC" >> str_$i.xsf
cat primvec.temp >> str_$i.xsf
echo "PRIMCOORD" >> str_$i.xsf
echo "$num_atom 1" >> str_$i.xsf
end=`echo "$i*$num_atom" |bc `
head -n $end pos.temp|tail -n $num_atom >pos_i.temp
paste type.temp pos_i.temp >> str_$i.xsf
mv str_$i.xsf $out-$i.xsf
done
rm *.temp
if [ ! -d "struc" ]; then
mkdir struc
fi
mv *xsf struc
Convert POSCAR to lammps input files¶
#!/bin/awk -f
#############################################################################
# #
# GNU License - Author: Pavlin D. Mitev 2012-10-09 #
# Version 0.1 #
# #
# Converts VASP POSCAR to LAMMPS imput structure. #
# It can handle non-orthogonal simulation boxes. #
# Reads both v4.6 and 5.2 (atom labels) POSCAR. #
# Accepts "Direct" and "Cartesian" coordinates in the POSCAR. #
# #
# Syntax: #
# VASP-poscar2lammps.awk POSCAR > structure.lammps #
# #
#############################################################################
BEGIN{
pi=3.14159265358979;rad2deg=180./pi;
if(ARGC<=1) { print "Syntax: \n VASP-poscar2res.awk POSCAR ..."; ex=1;exit}
for(i=2;i<=ARGC;i++){typeT[i-1]=ARGV[i];}
ARGC=2;
# Read header, scale and the basis -----------
getline; title=$0
getline; scale=$1
getline; h1[1]=$1*scale; h1[2]=$2*scale; h1[3]=$3*scale;
getline; h2[1]=$1*scale; h2[2]=$2*scale; h2[3]=$3*scale;
getline; h3[1]=$1*scale; h3[2]=$2*scale; h3[3]=$3*scale;
a=norm(h1); b=norm(h2); c=norm(h3); # Length of the basis vectors
alpha= angle(h2,h3); beta= angle(h1,h3); gamma= angle(h1,h2); # Angles in degree
alphar= alpha/rad2deg; betar= beta/rad2deg; gammar= gamma/rad2deg; # Angles in radians
# Check for labels -------------------------
getline;
if ($1*1 != $1) {
for(i=1;i<=NF;i++){typeT[i]=$i;}
getline;
}
for(i=1;i<=NF;i++) {type[i]=$i; natoms=natoms+$i} ntypes=NF;
# Advance to the coordinate section
while((tolower($0) !~ "direct")&&(tolower($0) !~ "cart")) getline;
if (tolower($0) ~ "direct") fractional=1; # Fractional format identified
# Rotation of the matrix to comply with LAMMPS standards ========================
p_a= sqrt(h1[1]**2 + h1[2]**2 + h1[3]**2);
p_b= sqrt(h2[1]**2 + h2[2]**2 + h2[3]**2);
p_c= sqrt(h3[1]**2 + h3[2]**2 + h3[3]**2);
lx= p_a;
p_xy= p_b * cos(gammar);
p_xz= p_c * cos(betar);
ly= sqrt(p_b**2 - p_xy**2);
p_yz= (p_b*p_c*cos(alphar)-p_xy*p_xz)/(ly);
lz= sqrt(p_c**2 - p_xz**2 - p_yz**2);
# The new basis H matrix ------------------------
H1[1]=lx; H1[2]= 0.000; H1[3]= 0.00;
H2[1]=p_xy; H2[2]= ly; H2[3]= 0.00;
H3[1]=p_xz; H3[2]= p_yz; H3[3]= lz;
# Matrix for conversion from cartesian to fractional in the old basis set (if necessary) ----------
cfv= sqrt(1. -cos(alphar)**2 -cos(betar)**2 -cos(gammar)**2 + 2.*cos(alphar)*cos(betar)*cos(gammar));
cf1[1]= 1./a; cf1[2]= -cos(gammar)/(a*sin(gammar)); cf1[3]= (cos(alphar)*cos(gammar)-cos(betar))/(a*cfv*sin(gammar));
cf2[1]= 0.00; cf2[2]= 1./(b*sin(gammar)); cf2[3]= (cos(betar)*cos(gammar)-cos(alphar)i)/(b*cfv*sin(gammar));
cf3[1]= 0.00; cf3[2]= 0.00; cf3[3]= sin(gammar)/(c*cfv);
# ===============================================================================
print "# Converted from POSCAR to lammps format"
print ""
print natoms" atoms"
print ntypes" atom types"
print ""
printf "0.000000 %10.6f xlo xhi\n",lx
printf "0.000000 %10.6f ylo yhi\n",ly
printf "0.000000 %10.6f zlo zhi\n",lz
print ""
printf "%10.6f %10.6f %10.6f xy xz yz\n", p_xy,p_xz,p_yz
print ""
print "Atoms"
print ""
iatom=0;
for(k=1;k<=ntypes;k++){
for(i=1;i<=type[k];i++){
getline
iatom++
x=$1; y=$2; z=$3;
if (!fractional){
xx=x*cf1[1]+y*cf1[2]+z*cf1[3];
yy=x*cf2[1]+y*cf2[2]+z*cf2[3];
zz=x*cf3[1]+y*cf3[2]+z*cf3[3];
x= xx; y= yy; z= zz;
}
xx=x*H1[1]+y*H2[1]+z*H3[1];
yy=x*H1[2]+y*H2[2]+z*H3[2];
zz=x*H1[3]+y*H2[3]+z*H3[3];
printf"%4i %-4s 0 %7f %7f %7f\n",iatom, k, xx,yy,zz
}
}
}
function asin(a) { return atan2(a,sqrt(1-a*a)) }
function acos(a) { return pi/2-asin(a) }
function norm(x) {return (sqrt(x[1]*x[1]+x[2]*x[2]+x[3]*x[3]));}
function dotprod (x,y) {return ( x[1]*y[1] + x[2]*y[2] + x[3]*y[3] );}
function angle (v1,v2) {
myacos = dotprod(v1,v2)/norm(v1)/norm(v2);
if (myacos>1.0) myacos = 1.0;
if (myacos<-1.0) myacos = -1.0;
return(acos(myacos)*180.0/3.14159265358979);
}
pos2lmp.awk POSCAR >data
awk '{if (NR!=9){
if(NR==11) {
printf("\nMasses\n\n1 64.0000\n2 12.01115\n\n");}
else print $0;
}
}' data >data.Cu
rm data
Convert POSCAR to XSF files¶
#lipai@mail.ustc.edu.cn
filename=$1
pos2xyz.pl $filename
echo "CRYSTAL" >> $filename.xsf
echo "PRIMVEC" >> $filename.xsf
head -5 $filename|tail -3 >> $filename.xsf
echo "PRIMCOORD" >> $filename.xsf
num=`head -1 $filename.xyz`
echo $num 1 >> $filename.xsf
awk '{if(NR>2) print $0}' $filename.xyz >> $filename.xsf
rm $filename.xyz
Calculate species numbers in a large configuration¶
This python script calculate number of speces from dump.trajlammps file
#lipai@mail.ustc.edu.cn
import numpy as np
class atom():
"""
Class for representing an atom
Parameters:
index: index of this atom
pos: position of this atom
sp_index: index of its species
atomnum: number of carbon atoms in its species
atoms_list: index of all species
if index!=sp_index, the info of atomnum and atoms_list is useless
Method:
join(index_j): add j atom into this species
"""
def __init__(self,index,pos):
self.index=index
self.pos=pos
self.sp_index=index
self.atomnum=1
self.atoms_list=[]
self.atoms_list.append(index)
def extend_sp(self,index_j):
global atoms
if self.index != self.sp_index:
print("error")
raise RunError
at=atoms[index_j]
self.atomnum=self.atomnum+at.atomnum
self.atoms_list.extend(at.atoms_list)
at.update(self.sp_index)
def update(self,index_j):
global atoms
for i in self.atoms_list:
atoms[i].sp_index=index_j;
def dist(i,j):
return np.linalg.norm(i.pos-j.pos)
def sp_join(i,j):
global atoms
if i.sp_index==j.sp_index:
return
if i.sp_index < j.sp_index:
atoms[i.sp_index].extend_sp(j.sp_index)
else:
atoms[j.sp_index].extend_sp(i.sp_index)
def calculate(pos):
global atoms
atoms=[]
for i in range(0,pos.shape[0]):
atoms.append(atom(i,pos[i]))
for i in range(0,len(atoms)):
for j in range(i,len(atoms)):
if dist(atoms[i],atoms[j])<1.7:
sp_join(atoms[i],atoms[j])
num=[0,0,0,0,0]
for i in atoms:
if i.index==i.sp_index:
if i.atomnum<5:
num[i.atomnum-1]+=1
return num
def skip(traj):
if len(traj.readline())==0:
return 0
for i in range(8):
traj.readline()
return 1
def get_atom_num(traj):
for i in range(3):
traj.readline()
num_atom=int(traj.readline().rstrip('\n'))
return num_atom
def write_info(sp):
for filename in filenames:
if len(filename)>0:
sp.write(str(num[0])+'\t'+str(num[1])+'\t'+str(num[2])+'\t'+str(num[3])+'\t'+str(num[4])+'\n')
atoms=[]
fin=open('dump.lammpstrj','r')
fout=open('species','w')
atom_num=get_atom_num(fin)
fin.seek(0)
while skip(fin):
pos=[]
for i in range(atom_num):
cont=fin.readline().rstrip('\n').split()
if len(cont) !=8:
print("len(cout): "+str(len(cont)))
print("cout.{}".format(cont))
print(str(len(pos)))
exit()
else:
if cont[1]=="2":
pos.append(cont)
pos=np.array(pos)
pos=pos[:,2:5]
pos=pos.astype(np.float64)
num=calculate(pos)
fout.write(str(num[0])+'\t'+str(num[1])+'\t'+str(num[2])+'\t'+str(num[3])+'\t'+str(num[4])+'\n')
fin.close()
fout.close()
Extract poscar files from traj.xyz¶
This python script extract poscar files from traj.xyz lattice constants are getton from CuC.xsf file
#lipai@mail.ustc.edu.cn
import os
f0=open("CuC.xsf")
f0.readline()
f0.readline()
x=f0.readline()
y=f0.readline()
z=f0.readline()
f0.close()
num_Cu=72
f1=open("traj.xyz")
ordi=1
while 1:
a=f1.readline()
if not a:
break
num=int(a.split()[0])
fileout="POSCAR-"+str(ordi)
f2=open(fileout,"w")
f2.write("Cu & C system from GCMC simulation\n")
f2.write("1.0\n")
f2.write(x)
f2.write(y)
f2.write(z)
f2.write("Cu C\n")
f2.write(str(num_Cu)+" "+str(num-num_Cu))
f2.write("\nCart\n")
ordi+=1
f1.readline()
for i in range(0,num):
a=f1.readline().split()
f2.write(a[1]+"\t"+a[2]+"\t"+a[3]+"\n")
f2.close()
f1.close()
os.system('mkdir poscar')
os.system('mv POSCAR* poscar')
Python Scripts¶
This page contains python scripts
Scripts¶
Calculate ionic conductivity in solid¶
This script calculate diffusivity and conductivity from vasprun.xml file
import matplotlib.pyplot as plt
#import json
import collections
from pymatgen import Structure
from pymatgen.analysis.diffusion_analyzer import DiffusionAnalyzer, \
get_arrhenius_plot, get_extrapolated_conductivity
from pymatgen_diffusion.aimd.pathway import ProbabilityDensityAnalysis
from pymatgen_diffusion.aimd.van_hove import VanHoveAnalysis
temperatures = [500, 1000, 1500]
#analyzers = collections.OrderedDict()
#for temp in temperatures:
# with open("%d.json" % temp) as f:
# d = json.load(f)
# analyzers[temp] = DiffusionAnalyzer.from_dict(d)
#
#plt = analyzers[1000].get_msd_plot()
#title = plt.title("1000K", fontsize=24)
#
#diffusivities = [d.diffusivity for d in analyzers.values()]
#plt = get_arrhenius_plot(temperatures, diffusivities)
#
#rts = get_extrapolated_conductivity(temperatures, diffusivities,
# new_temp=300, structure=analyzers[800].structure,
# species="Li")
#print("The Li ionic conductivity for Li6PS5Cl at 300 K is %.4f mS/cm" % rts)
a = ["500/vasprun.xml"]
b = ["1000/vasprun.xml"]
c = ["1500/vasprun.xml"]
analyzer_500 = DiffusionAnalyzer.from_files(a, specie="Li", smoothed=False)
analyzer_1000 = DiffusionAnalyzer.from_files(b, specie="Li", smoothed=False)
analyzer_1500 = DiffusionAnalyzer.from_files(c, specie="Li", smoothed=False)
diffusivities=[analyzer_500.diffusivity,analyzer_1000.diffusivity,analyzer_1500.diffusivity]
plt = get_arrhenius_plot(temperatures, diffusivities)
plt.savefig("Li_Arrhenius_Plot.png")
rts = get_extrapolated_conductivity(temperatures, diffusivities,
new_temp=300, structure=analyzer_500.structure,
species="Li")
print("The Li ionic conductivity for anti-spinel Li3OBr at 300 K is %.4f mS/cm" % rts)
Get MSD (mean square displacement) from MD trajectory¶
The first script convert XDATCAR to XYZ file The second script extract MSD info from XYZ file
import numpy as np
franum=0
count=0
latt=[]
axyz=[]
ele=[]
etnum=[]
#fixflag=input(" Seletive Dynamics is T or F: ")
fixflag=False
tlnum=0
pf=open("XDATCAR","r")
fileend=0
tmpstr=pf.readline()
count=count+1
tmpstr=pf.readline()
count=count+1
#### Begin Read Lattice Vector #####
latt.append([])
for j in range(0,3):
tl=pf.readline()
count=count+1
tlst=tl.split()
tlf=[float(x) for x in tlst[0:3]]
latt[franum].append(tlf)
print latt
#### Begin Read Element Num #####
tmpstr=pf.readline()
count=count+1
[ele.append(x) for x in tmpstr.split()]
print ele
tnumstr=pf.readline()
print "read ele num"
print tnumstr
count=count+1
etnum=[int(x) for x in tnumstr.split()]
print etnum
totnum=sum(etnum)
if fixflag:
tmpstr=pf.readline()
count=count+1
#### Begin Read Atom Fraction Coordination #####
while not fileend:
tmpstr=pf.readline()
if tmpstr == '':
print "file read done"
fileend=1
break
count=count+1
axyz.append([])
for j in range(0,totnum):
txyz=pf.readline()
if txyz == '':
print "Error:Missing Atom In Frame %d"%franum
fileend=1
break
count=count+1
taxyz=txyz.split()
xyz=[float(x) for x in taxyz[0:3]]
axyz[franum].append(xyz)
franum=franum+1
LATTICE=np.array(latt)
ATXYZ=np.array(axyz)
print ATXYZ[0]
CXYZ=[]
outfile="LiXY"
of=open(outfile,"w")
for i in range(0,franum,10):
# ATXYZ=np.array(axyz[i])
# LATT=LATTICE[i]
# tCXYZ=np.dot(ATXYZ,LATT)
# CXYZ.append(tCXYZ)
# print >> of," %d "%(enum[i][2])
# print >> of," %s %d %s "%("### ",i+1," ####")
for j in range(0,etnum[0]):
print >> of," %f %f "%(ATXYZ[i][j][0],ATXYZ[i][j][1])
of.close()
# The MIT License (MIT)
#
# Copyright (c) 2014 Muratahan Aykol
#
# Permission is hereby granted, free of charge, to any person obtaining a copy
# of this software and associated documentation files (the "Software"), to deal
# in the Software without restriction, including without limitation the rights
# to use, copy, modify, merge, publish, distribute, sublicense, and/or sell
# copies of the Software, and to permit persons to whom the Software is
# furnished to do so, subject to the following conditions:
#
# The above copyright notice and this permission notice shall be included in
# all copies or substantial portions of the Software.
#
# THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR
# IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
# FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE
# AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
# LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM,
# OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN
# THE SOFTWARE
import numpy as np
from copy import deepcopy
# This function reads an XYZ file and a list of lattice vectors L = [x,y,z] and gives MSD + unwrapped coordinates
def MSD(xyz_file,L):
a = []; l = [];
a.append(L[0]); a.append(L[1]); a.append(L[2]); #basis vectors in cartesian coords
l.append(np.sqrt(np.dot(a[0],a[0]))); l.append(np.sqrt(np.dot(a[1],a[1]))); l.append(np.sqrt(np.dot(a[2],a[2]))); #basis vector lengths
file = open(xyz_file, 'r')
recorder = open("msd.out", 'w')
coord_rec = open("unwrapped.xyz", 'w')
origin_list = [] # Stores the origin as [element,[coords]]
prev_list = [] # Stores the wrapped previous step
unwrapped_list = [] # Stores the instantenous unwrapped
msd = [] #Stores atom-wise MSD Stores msd as [msd]
msd_dict ={} #Stores element-wise MSD
msd_lattice = []
msd_dict_lattice ={}
element_list = [] # element list
element_dict = {} # number of elements stored
content = file.readline()
N = int(content)
for i in range(N):
msd.append(np.float64('0.0'))
msd_lattice.append([0.0, 0.0, 0.0 ])
file.readline()
step = 0
while True:
step += 1
# Get and store the origin coordinates in origin_dict at first step
if step == 1:
for i in xrange(N):
t = file.readline().rstrip('\n').split()
element = t[0]
if element not in element_list:
element_list.append(element)
if element not in element_dict:
element_dict[element] = 1.0
else:
element_dict[element] += 1.0
coords = np.array( [ float(s) for s in t[1:] ] )
origin_list.append([element,coords])
# Copy the first set of coordinates as prev_dict and unwrapped
unwrapped_list = deepcopy(origin_list)
prev_list = deepcopy(origin_list)
recorder.write("step ")
for element in element_list:
recorder.write(element+" ")
recorder.write("\n")
# Read wrapped coordinates into wrapped_dict
content = file.readline()
if len(content) == 0:
print "\n---End of file---\n"
break
N = int(content)
file.readline()
wrapped_list = [] # Erease the previous set of coordinates
for i in xrange(N):
t = file.readline().rstrip('\n').split()
element = t[0]
coords = np.array( [ float(s) for s in t[1:] ] )
wrapped_list.append([element,coords])
coord_rec.write(str(N)+ "\ncomment\n")
# Unwrap coodinates and get MSD
for atom in range(N):
msd[atom] = 0.0
coord_rec.write(wrapped_list[atom][0])
# decompose wrapped atom coordinates to onto lattice vectors:
w1 = wrapped_list[atom][1][0]
w2 = wrapped_list[atom][1][1]
w3 = wrapped_list[atom][1][2]
# decompose prev atom coordinates to onto lattice vectors:
p1 = prev_list[atom][1][0]
p2 = prev_list[atom][1][1]
p3 = prev_list[atom][1][2]
#get distance between periodic images and use the smallest one
if np.fabs(w1 - p1) > 0.5:
u1 = w1 - p1 - np.sign(w1 - p1)
else:
u1 = w1 - p1
if np.fabs(w2 - p2) > 0.5:
u2 = w2 - p2 - np.sign(w2 - p2)
else:
u2 = w2 - p2
if np.fabs(w3 - p3) > 0.5:
u3 = w3 - p3 - np.sign(w3 - p3)
else:
u3 = w3 - p3
#add unwrapped displacements to unwrapped coords
unwrapped_list[atom][1][0] += u1
unwrapped_list[atom][1][1] += u2
unwrapped_list[atom][1][2] += u3
uw = unwrapped_list[atom][1][0]*a[0] + unwrapped_list[atom][1][1]*a[1] +unwrapped_list[atom][1][2]*a[2]
ol = origin_list[atom][1][0]*a[0] + origin_list[atom][1][1]*a[1] + origin_list[atom][1][2]*a[2]
msd[atom] = np.linalg.norm(uw-ol)**2
msd_lattice[atom] = [np.linalg.norm(uw[0]-ol[0])**2,np.linalg.norm(uw[1]-ol[1])**2,np.linalg.norm(uw[2]-ol[2])**2]
coord_rec.write(" " + np.array_str(uw).replace("[","").replace("]",""))
coord_rec.write("\n")
prev_list = [] # Store current wrapped coordinates for the next step
prev_list = deepcopy(wrapped_list)
# record msd
recorder.write(str(step) + " ")
for el in element_list:
msd_dict[el] = 0.0
msd_dict_lattice[el]=[0.,0.,0.]
for atom in range(len(msd)):
msd_dict[wrapped_list[atom][0]] += msd[atom]/element_dict[wrapped_list[atom][0]]
for i in range(3):
msd_dict_lattice[wrapped_list[atom][0]][i] += msd_lattice[atom][i]/element_dict[wrapped_list[atom][0]]
for el in element_list:
recorder.write(str(msd_dict[el])+ " " + str(msd_dict_lattice[el][0])+ " " + str(msd_dict_lattice[el][1])+ " " + str(msd_dict_lattice[el][2])+ " ")
recorder.write("\n")
if step % 10 == 0:
print step
recorder.close()
file.close()
coord_rec.close()
def read_lat_vec():
lat_file = open('lattice.vectors','r')
line = []
for i in xrange(3):
line.append([float(x) for x in lat_file.readline().rstrip('\n').split()])
print line[i]
lattice = np.array([line[0],line[1],line[2]])
return lattice
#You can read the lattice vectors from lattice.vector file
lattice = read_lat_vec()
#Or define the lattice vector manually as in
#lattice =np.array([[-12.181156,-4.306689,7.459404],[0.000000,-12.920067,7.459404],[0.000000,0.000000,14.918808]])
#Run the MSD calculator with XDATCAR_fract.xyz and the lattice vector defined above
MSD("XDATCAR_fract.xyz",lattice)
Plot Phasediagram¶
This python script use ase to plot phasediagram
#lipai@mail.ustc.edu.cn
from ase.phasediagram import PhaseDiagram
u=[0.0,-2.92,-3.10,-3.20]
Li=0
Br=0.0
O2=0.0
Br2O3=-0.209*5
for u_Li in u:
Li2O=-2.071*3-2*u_Li
Li2O2=-1.65*4-2*u_Li
LiBr=-1.573*2-u_Li
#Li3OBr=Li2O+LiBr+0.075*5
Li32O11Br10=11*Li2O+10*LiBr+0.059*53
refs = [('Br',Br),
('O2',O2),
('Li',u_Li),
('Br2O3',Br2O3),
('Li2O',Li2O),
('Li2O2',Li2O2),
('LiBr',LiBr),
# ('Li3OBr',Li3OBr),
('Li32O11Br10',Li32O11Br10)]
pd = PhaseDiagram(refs)
fig=pd.plot(show=False)
fig.savefig("try"+str(u_Li)+".png")
#pd.savefig(str(u_Li)+".png")
import matplotlib.pyplot as plt
u=[0.0,-2.92,-3.10,-3.20]
for u_Li in u:
O2=0.0
Li2O=(-2.071*3-2*u_Li)/1
Li2O2=(-1.65*4-2*u_Li)/2
Br2O3=(-0.209*5)/5
Br=0.0
LiBr=(-1.573*2-u_Li)/1
Li32O11Br10=(11*Li2O+10*LiBr+0.059*53)/21
x_Br=[1,0,2/5.0,0,0,1,10/21.0]
y=[Br,O2,Br2O3,Li2O,Li2O2,LiBr,Li32O11Br10]
name=["Br","O2","Br2O3","Li2O","Li2O2","LiBr","Li32O11Br10"]
plt.plot(x_Br,y,'sr')
for i,j,k in zip(x_Br,y,name):
#cname=
if k=="Li2O" :
plt.text(i-0.09,j+0.01,k,color='g')
else:
plt.text(i+0.02,j+0.01,k,color='r')
plt.xlim([-0.1,1.1])
plt.axis("off")
plt.savefig(str(u_Li)+".png")
plt.clf()
print "x_Br y"
for i,j in zip(x_Br,y):
print i,j
Plot plain figures¶
Usage: python plot.py [filename] column_x column_y1 [column_y2 …]
plot.py for plot data and save fig¶
#lipai@mail.ustc.edu.cn
import sys
import matplotlib as mpl
import matplotlib.pyplot as plt
import numpy as np
mpl.rcParams['agg.path.chunksize']=10000
arg_num=len(sys.argv)
if arg_num==1:
print("no inputs!!!\n")
exit()
filename=sys.argv[1]
data=np.loadtxt(filename,skiprows=1)
colors=['k','r','g','b','c','m']
labels=['DFT','NN']
plt.figure(figsize=(10,4))
if arg_num==2:
plt.plot(data[:,0],data[:,1])
elif arg_num>3:
yn=arg_num-3
if int(sys.argv[2]) == 0:
# x=np.arange(data.shape)
# else:
# x=np.arange(data.shape[0])
x=np.arange(list(data.shape)[0])
for i in range(yn):
if i==0:
alpha=1
else:
alpha=0.7
if len(data.shape)==1:
plt.plot(x,data[:],color=colors[i],linewidth=2,label=labels[i],alpha=alpha)
else:
plt.plot(x,data[:,int(sys.argv[i+3])-1],color=colors[i],linewidth=2,label=labels[i],alpha=alpha)
else:
for i in range(yn):
plt.plot(data[:,int(sys.argv[2])-1],data[:,int(sys.argv[i+3])-1],color=colors[i],linewidth=2,label=labels[i])
else:
print("wrong input arguments")
plt.ylabel('Energy/eV')
#plt.ylim(0,5)
#plt.axhline(0)
plt.legend(loc=0)
plt.savefig("temp.jpg",dpi=300)
#plt.savefig("temp.jpg")
#plt.show()
Windows:¶
断网自动重连¶
1.右键我的电脑—管理 2.系统工具–任务计划程序–创建基本任务向导 3.名称：断网链接；描述：断网重连；点下一步 4.触发器：选“当特定事件被记录时”；下一步 5.日志：系统；源：rasman；事件：20268；下一步 6.直接下一步 7.程序或脚本：rasdial 宽带连接 用户名 密码
Tutorials:¶
Vim skills¶
代码折叠¶
vim 提供 6中折叠方式
manual 手工定义折叠
indent 更多的缩进表示更高级别的折叠
expr 用表达式来定义折叠
syntax 用语法高亮来定义折叠
diff 对没有更改的文本进行折叠
marker 对文中的标志折叠
可用选项 ‘foldmethod’ 来设定折叠方式：set fdm=syntax
- 折叠打开与折合
zc 折叠
zC 对所在范围内所有嵌套的折叠点进行折叠
zo 展开折叠
zO 对所在范围内所有嵌套的折叠点展开
[z 到当前打开的折叠的开始处。
]z 到当前打开的折叠的末尾处。
zj 向下移动。到达下一个折叠的开始处。关闭的折叠也被计入。
zk 向上移动到前一折叠的结束处。关闭的折叠也被计入。
My .vimrc file¶
"set nobackup " 关闭自动备份功能，backup自动备份
set scrolloff=5 " 光标移动到buffer的顶部和底部时保持*行距离
set nu
colorscheme desert
syntax enable
syntax on
set nowrap
set expandtab
set shiftwidth=4
set tabstop=4
let fortran_fold=1
set foldmethod=syntax
set foldlevelstart=10
set mouse=a
set paste
A ReStructuredText Primer¶
The text below contains links that look like “(quickref)”. These are relative links that point to the Quick reStructuredText user reference. If these links don’t work, please refer to the master quick reference document.
Note
This document is an informal introduction to reStructuredText. The `What Next?`_ section below has links to further resources, including a formal reference.
Structure¶
From the outset, let me say that “Structured Text” is probably a bit of a misnomer. It’s more like “Relaxed Text” that uses certain consistent patterns. These patterns are interpreted by a HTML converter to produce “Very Structured Text” that can be used by a web browser.
The most basic pattern recognised is a paragraph (quickref). That’s a chunk of text that is separated by blank lines (one is enough). Paragraphs must have the same indentation – that is, line up at their left edge. Paragraphs that start indented will result in indented quote paragraphs. For example:
This is a paragraph. It's quite
short.
This paragraph will result in an indented block of
text, typically used for quoting other text.
This is another one.
Results in:
This is a paragraph. It’s quite short.
This paragraph will result in an indented block of text, typically used for quoting other text.This is another one.
Text styles¶
(quickref)
Inside paragraphs and other bodies of text, you may additionally mark
text for italics with “*italics*” or bold with
“**bold**”. This is called “inline markup”.
If you want something to appear as a fixed-space literal, use
“``double back-quotes``”. Note that no further fiddling is done
inside the double back-quotes – so asterisks “*” etc. are left
alone.
If you find that you want to use one of the “special” characters in
text, it will generally be OK – reStructuredText is pretty smart.
For example, this lone asterisk * is handled just fine, as is the
asterisk in this equation: 5*6=30. If you actually
want text *surrounded by asterisks* to not be italicised, then
you need to indicate that the asterisk is not special. You do this by
placing a backslash just before it, like so “\*” (quickref), or
by enclosing it in double back-quotes (inline literals), like this:
``*``
Tip
Think of inline markup as a form of (parentheses) and use it the same way: immediately before and after the text being marked up. Inline markup by itself (surrounded by whitespace) or in the middle of a word won’t be recognized. See the markup spec for full details.
Lists¶
Lists of items come in three main flavours: enumerated, bulleted and definitions. In all list cases, you may have as many paragraphs, sublists, etc. as you want, as long as the left-hand side of the paragraph or whatever aligns with the first line of text in the list item.
Lists must always start a new paragraph – that is, they must appear after a blank line.
- enumerated lists (numbers, letters or roman numerals; quickref)
Start a line off with a number or letter followed by a period “.”, right bracket “)” or surrounded by brackets “( )” – whatever you’re comfortable with. All of the following forms are recognised:
1. numbers A. upper-case letters and it goes over many lines with two paragraphs and all! a. lower-case letters 3. with a sub-list starting at a different number 4. make sure the numbers are in the correct sequence though! I. upper-case roman numerals i. lower-case roman numerals (1) numbers again 1) and again
Results in (note: the different enumerated list styles are not always supported by every web browser, so you may not get the full effect here):
- numbers
upper-case letters and it goes over many lines
with two paragraphs and all!
- lower-case letters
- with a sub-list starting at a different number
- make sure the numbers are in the correct sequence though!
- upper-case roman numerals
- lower-case roman numerals
- numbers again
- and again
- bulleted lists (quickref)
Just like enumerated lists, start the line off with a bullet point character - either “-“, “+” or “*”:
* a bullet point using "*" - a sub-list using "-" + yet another sub-list - another item
Results in:
- a bullet point using “*”
- a sub-list using “-“
- yet another sub-list
- another item
- a sub-list using “-“
- a bullet point using “*”
- definition lists (quickref)
Unlike the other two, the definition lists consist of a term, and the definition of that term. The format of a definition list is:
what Definition lists associate a term with a definition. *how* The term is a one-line phrase, and the definition is one or more paragraphs or body elements, indented relative to the term. Blank lines are not allowed between term and definition.
Results in:
- what
- Definition lists associate a term with a definition.
- how
- The term is a one-line phrase, and the definition is one or more paragraphs or body elements, indented relative to the term. Blank lines are not allowed between term and definition.
Preformatting (code samples)¶
(quickref)
To just include a chunk of preformatted, never-to-be-fiddled-with
text, finish the prior paragraph with “::”. The preformatted
block is finished when the text falls back to the same indentation
level as a paragraph prior to the preformatted block. For example:
An example::
Whitespace, newlines, blank lines, and all kinds of markup
(like *this* or \this) is preserved by literal blocks.
Lookie here, I've dropped an indentation level
(but not far enough)
no more example
Results in:
An example:
Whitespace, newlines, blank lines, and all kinds of markup (like *this* or \this) is preserved by literal blocks. Lookie here, I've dropped an indentation level (but not far enough)no more example
Note that if a paragraph consists only of “::”, then it’s removed
from the output:
::
This is preformatted text, and the
last "::" paragraph is removed
Results in:
This is preformatted text, and the
last "::" paragraph is removed
| Author: | Richard Jones |
|---|
Links:¶
Galleries:¶
Supports:¶
This website is created using Sphinx