vasp_tools/create_lattice.c

   1 /*
   2  * create lattice (for usage in vasp POSCAR file)
   3  *
   4  * author: Frank Zirkelbach <frank.zirkelbach@physik.uni-augsburg.de>
   5  *
   6  */
   7
   8 #include <stdio.h>
   9 #include <stdlib.h>
  10 #include <string.h>
  11 #include <math.h>
  12 #include "../math/math.h"
  13
  14 /*
  15  * lattice types:
  16  *
  17  * there is fcc/diamond/zincblende only!
  18  *
  19  * basis:
  20  *
  21  * 0: <0.5 0.5 0> <0 0.5 0.5> <0.5 0 0.5> contains 2 atoms
  22  * 1: <0.5 -0.5 0> <0.5 0.5 0> <0 0 1> contains 4 atoms
  23  * 2: <1 0 0> <0 1 0> <0 0 1> contains 8 atoms
  24  *
  25  * atoms:
  26  *
  27  * - type 0:
  28  *   1st fcc: <0 0 0>
  29  *   2nd fcc: <0.25 0.25 0.25>
  30  *   -> diag: <0.25 0.25 0.25>
  31  *
  32  * - type 1:
  33  *   1st fcc: <0 0 0> <0.5 0.5 0.5>
  34  *   2nd fcc: <0.5 0 0.25> <1.0 or 0.0 0.5 0.75>
  35      -> diag: <0.5 0 0.25>
  36  *
  37  * - type 2:
  38  *   1st fcc: <0 0 0> <0.5 0.5 0> <0.5 0 0.5> <0 0.5 0.5>
  39  *   2nd fcc: <0.25 0.25 0.25> <0.75 0.75 0.25> <0.75 0.25 0.75>
  40  *            <0.25 0.75 0.75>
  41  *   -> diag: <0.25 0.25 0.25>
  42  *
  43  */
  44
  45 int main(int argc,char **argv) {
  46
  47         int i,j,k,l,cnt,estimated;
  48         int x,y,z;
  49         t_3dvec basis[3];
  50         t_3dvec o[3];
  51         t_3dvec dia;
  52         t_3dvec r,h;
  53         char type,fccdia;
  54
  55         if(argc!=6) {
  56                 printf("usage: %s type fcc/dia lx ly lz\n",argv[0]);
  57                 printf("  basis types: 0, 1, 2 (see code)\n");
  58                 printf("  fcc/dia: 0=fcc, 1=dia\n");
  59                 return -1;
  60         }
  61
  62         type=argv[1][0];
  63         fccdia=argv[2][0];
  64
  65         x=atoi(argv[3]);
  66         y=atoi(argv[4]);
  67         z=atoi(argv[5]);
  68
  69         dia.x=0.25;
  70         dia.y=0.25;
  71         dia.z=0.25;
  72
  73         if(type=='1') {
  74                 dia.x=0.5;
  75                 dia.y=0.0;
  76                 dia.z=0.25;
  77         }
  78
  79         memset(basis,0,3*sizeof(t_3dvec));
  80         memset(o,0,3*sizeof(t_3dvec));
  81
  82         if(type=='0') {
  83                 basis[0].x=0.5;
  84                 basis[0].y=0.5;
  85                 basis[1].y=0.5;
  86                 basis[1].z=0.5;
  87                 basis[2].x=0.5;
  88                 basis[2].z=0.5;
  89         }
  90
  91         if(type=='1') {
  92                 basis[0].x=0.5;
  93                 basis[0].y=-0.5;
  94                 basis[1].x=0.5;
  95                 basis[1].y=0.5;
  96                 basis[2].z=1.0;
  97                 o[0].x=0.5;
  98                 o[0].y=0.5;
  99                 o[0].z=0.5;
 100         }
 101
 102         if(type=='2') {
 103                 basis[0].x=1.0;
 104                 basis[1].y=1.0;
 105                 basis[2].z=1.0;
 106                 o[0].x=0.5;
 107                 o[0].y=0.5;
 108                 o[1].x=0.5;
 109                 o[1].z=0.5;
 110                 o[2].y=0.5;
 111                 o[2].z=0.5;
 112         }
 113
 114         cnt=0;
 115
 116         estimated=1;
 117         if(fccdia=='1') estimated*=2;
 118         if(type=='1') estimated*=2;
 119         if(type=='2') estimated*=4;
 120         estimated*=x*y*z;
 121
 122         // print POSCAR 'header'
 123
 124         printf("cubic diamond\n");
 125         printf(" 5.429\n");
 126
 127         v3_scale(&h,&basis[0],x);
 128         printf(" %.5f %.5f %.5f\n",h.x,h.y,h.z);
 129         v3_scale(&h,&basis[1],y);
 130         printf(" %.5f %.5f %.5f\n",h.x,h.y,h.z);
 131         v3_scale(&h,&basis[2],z);
 132         printf(" %.5f %.5f %.5f\n",h.x,h.y,h.z);
 133
 134         printf(" %d\n",estimated);
 135         printf("selective dynamics\n");
 136         printf("direct\n");
 137
 138         // now print the coordinates
 139
 140         for(i=0;i<x;i++) {
 141                 for(j=0;j<y;j++) {
 142                         for(k=0;k<z;k++) {
 143
 144         // first atom (all types)
 145         r.x=i;
 146         r.y=j;
 147         r.z=k;
 148
 149         printf(" %.5f %.5f %.5f T T T\n",r.x/x,r.y/y,r.z/z);
 150         cnt+=1;
 151
 152         // second atom (type 1)
 153         if(type=='1') {
 154                 v3_add(&h,&r,&(o[0]));
 155                 printf(" %.5f %.5f %.5f T T T\n",h.x/x,h.y/y,h.z/z);
 156                 cnt+=1;
 157         }
 158
 159         // face centered atoms (type 2)
 160         if(type=='2') {
 161                 for(l=0;l<3;l++) {
 162                         v3_add(&h,&r,&(o[l]));
 163                         printf(" %.5f %.5f %.5f T T T\n",h.x/x,h.y/y,h.z/z);
 164                         cnt+=1;
 165                 }
 166         }
 167
 168                         }
 169                 }
 170         }
 171
 172         // second fcc lattice
 173         if(fccdia=='1') {
 174
 175         for(i=0;i<x;i++) {
 176                 for(j=0;j<y;j++) {
 177                         for(k=0;k<z;k++) {
 178
 179         // first atom (all types)
 180         r.x=i;
 181         r.y=j;
 182         r.z=k;
 183
 184         // add the diagonal
 185         v3_add(&r,&r,&dia);
 186
 187         printf(" %.5f %.5f %.5f T T T\n",r.x/x,r.y/y,r.z/z);
 188         cnt+=1;
 189
 190         // second atom (type 1)
 191         if(type=='1') {
 192                 v3_add(&h,&r,&(o[0]));
 193                 printf(" %.5f %.5f %.5f T T T\n",h.x/x,h.y/y,h.z/z);
 194                 cnt+=1;
 195         }
 196
 197         // face centered atoms (type 2)
 198         if(type=='2') {
 199                 for(l=0;l<3;l++) {
 200                         v3_add(&h,&r,&(o[l]));
 201                         printf(" %.5f %.5f %.5f T T T\n",h.x/x,h.y/y,h.z/z);
 202                         cnt+=1;
 203                 }
 204         }
 205
 206                         }
 207                 }
 208         }
 209
 210
 211         }
 212
 213
 214         if(estimated!=cnt) {
 215                 printf("\nWARNING! counted %d ions and estimated %d ions\n\n",
 216                        cnt,estimated);
 217         }
 218
 219         return 0;
 220 }