added nlsop_client code
[physik/nlsop.git] / nlsop_client.c
1 /*
2  * client nlsop code
3  *
4  * author: frank zirkelbach (frank.zirkelbach@physik.uni-augsburg.de)
5  *
6  * this program tries helping to understand the amorphous depuration
7  * and recrystallization of SiCx while ion implantation at temperatures
8  * below 400 degree celsius.
9  * hopefully the program will simulate the stabilization of the
10  * selforganizing lamella structure in the observed behaviour.
11  *
12  * refs: 
13  *  - J. K. N. Lindner. Habil.Schrift, Universitaet Augsburg.
14  *  - Maik Haeberlen. Diplomarbeit, Universitaet Augsburg.
15  *
16  * Copyright (C) 2004 Frank Zirkelbach
17  *
18  * This program is free software; you can redistribute it and/or modify
19  * it under the terms of the GNU General Public License as published by
20  * the Free Software Foundation; either version 2 of the License, or
21  * (at your option) any later version.
22  *
23  * This program is distributed in the hope that it will be useful,
24  * but WITHOUT ANY WARRANTY; without even the implied warranty of
25  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
26  * GNU General Public License for more details.
27  *
28  * You should have received a copy of the GNU General Public License
29  * along with this program; if not, write to the Free Software
30  * Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA
31  *
32  */
33
34 #define _GNU_SOURCE
35 #include <stdio.h>
36 #include <stdlib.h>
37 #include <string.h>
38 #include <sys/types.h>
39 #include <sys/stat.h>
40 #include <fcntl.h>
41 #include <unistd.h>
42
43 #include "nlsop.h"
44 #include "dfbapi.h"
45 #include "random.h"
46
47 #include "../api/network.h"
48 #include "../api/event.h"
49 #include "../api/list.h"
50
51 #define MAKE_AMORPH(N) *(N)|=AMORPH
52 #define MAKE_CRYST(N) *(N)&=~AMORPH
53
54 int usage(char *prog)
55 {
56  puts("usage:");
57  printf("%s <server ip>\n",prog);
58  return 1;
59 }
60
61 /*
62  * nlsop internal functions
63  */
64
65 int sputter(d3_lattice *d3_l)
66 {
67  int i,size;
68  int offh,offl;
69
70  size=d3_l->max_x*d3_l->max_y;
71  offl=0;
72  offh=size;
73
74  for(i=0;i<d3_l->max_z-1;i++)
75  {
76   memcpy(d3_l->status+offl,d3_l->status+offh,size);
77   memcpy(d3_l->extra+offl,d3_l->extra+offh,size*sizeof(int));
78   offl=offh;
79   offh+=size;
80  }
81  memset(d3_l->status+offl,0,size);
82  memset(d3_l->extra+offl,0,size);
83
84  return 1;
85 }
86
87 int process_cell(d3_lattice *d3_l,u32 x,u32 y,u32 z,info *my_info,u32 nel_z)
88 {
89  unsigned char *thiz;
90  int *conc;
91  int i,j;
92  int off;
93  double p,q;
94
95  thiz=d3_l->status+x+y*d3_l->max_x+z*d3_l->max_x*d3_l->max_y;
96  conc=d3_l->extra+x+y*d3_l->max_x+z*d3_l->max_x*d3_l->max_y;
97  p=my_info->b*nel_z;
98  for(i=-(my_info->range);i<=my_info->range;i++)
99  {
100   for(j=-(my_info->range);j<=my_info->range;j++)
101   {
102    if(!(i==0 && j==0))
103    {
104     off=((x+d3_l->max_x+i)%d3_l->max_x)+((y+d3_l->max_y+j)%d3_l->max_x)*d3_l->max_x+z*d3_l->max_x*d3_l->max_y;
105     if(*(d3_l->status+off)&AMORPH) p+=my_info->s*(*(d3_l->extra+off))*URAND_MAX/(i*i+j*j);
106    } 
107   }
108  }
109  p+=*conc*my_info->c*URAND_MAX;
110  if(!(*thiz&AMORPH))
111  {
112   if(get_rand(URAND_MAX)<=p) MAKE_AMORPH(thiz);
113  } else
114  {
115   /* assume 1-p probability */
116   /* also look for neighbours ! */
117   q=(URAND_MAX-p)>0?URAND_MAX-p:0;
118   j=0;
119   j+=(*(d3_l->status+((x+d3_l->max_x+1)%d3_l->max_x)+y*d3_l->max_x+z*d3_l->max_x*d3_l->max_y)&AMORPH)?1:0;
120   j+=(*(d3_l->status+((x+d3_l->max_x-1)%d3_l->max_x)+y*d3_l->max_x+z*d3_l->max_x*d3_l->max_y)&AMORPH)?1:0;
121   j+=(*(d3_l->status+x+((y+1+d3_l->max_y)%d3_l->max_y)*d3_l->max_x+z*d3_l->max_x*d3_l->max_y)&AMORPH)?1:0;
122   j+=(*(d3_l->status+x+((y-1+d3_l->max_y)%d3_l->max_y)*d3_l->max_x+z*d3_l->max_x*d3_l->max_y)&AMORPH)?1:0;
123   j+=(*(d3_l->status+x+y*d3_l->max_x+((z+1+d3_l->max_z)%d3_l->max_z)*d3_l->max_x*d3_l->max_y)&AMORPH)?1:0;
124   j+=(*(d3_l->status+x+y*d3_l->max_x+((z-1+d3_l->max_z)%d3_l->max_z)*d3_l->max_x*d3_l->max_y)&AMORPH)?1:0;
125
126   p+=((q/6)*j);
127   if(get_rand(URAND_MAX)>p) MAKE_CRYST(thiz);
128  }
129  
130  return 1;
131 }
132
133 int distrib_c(d3_lattice *d3_l,info *my_info,int step,u32 rj_m,u32 *rj_g)
134 {
135  u32 x,y,z;
136  int i,j,k,c;
137  int offset,off;
138  int carry;
139
140  /* put one c ion somewhere in the lattice */
141  x=get_rand(d3_l->max_x);
142  y=get_rand(d3_l->max_y);
143  z=get_rand_reject(d3_l->max_z,rj_m,rj_g);
144  *(d3_l->extra+x+y*d3_l->max_x+z*d3_l->max_x*d3_l->max_y)+=1;
145  (my_info->cc)++;
146
147  if(step%my_info->diff_rate==0)
148  {
149
150  for(i=0;i<d3_l->max_x;i++)
151  {
152   for(j=0;j<d3_l->max_y;j++)
153   {
154    for(k=0;k<d3_l->max_z;k++)
155    {
156     offset=i+j*d3_l->max_x+k*d3_l->max_x*d3_l->max_y;
157     /* case amorph: amorph <- cryst diffusion */
158     if(*(d3_l->status+offset)&AMORPH)
159     {
160      for(c=-1;c<=1;c++)
161      {
162       if(c!=0)
163       {
164        off=((i+d3_l->max_x+c)%d3_l->max_x)+j*d3_l->max_x+k*d3_l->max_x*d3_l->max_y;
165        carry=0;
166        if(!(*(d3_l->status+off)&AMORPH)) carry=(int)(my_info->dr_ac*(*(d3_l->extra+off)));
167        if(carry!=0)
168        {
169         *(d3_l->extra+offset)+=carry;
170         *(d3_l->extra+off)-=carry;
171        }
172       }
173      }
174      for(c=-1;c<=1;c++)
175      {
176       if(c!=0)
177       {
178        off=i+((j+c+d3_l->max_y)%d3_l->max_y)*d3_l->max_x+k*d3_l->max_x*d3_l->max_y;
179        carry=0;
180        if(!(*(d3_l->status+off)&AMORPH)) carry=(int)(my_info->dr_ac*(*(d3_l->extra+off)));
181        if(carry!=0)
182        {
183         *(d3_l->extra+offset)+=carry; 
184         *(d3_l->extra+off)-=carry; 
185        }
186       }
187      }
188      if(my_info->z_diff)
189      {
190       if(k!=0)
191       {
192        off=i+j*d3_l->max_x+(k-1)*d3_l->max_x*d3_l->max_y;
193        carry=0;
194        if(!*(d3_l->status+off)&AMORPH) carry=(int)(my_info->dr_ac*(*(d3_l->extra+off)));
195        if(carry!=0)
196        {
197         *(d3_l->extra+off)-=carry;
198         *(d3_l->extra+offset)+=carry;
199        }
200       }
201       if(k!=d3_l->max_z-1)
202       {
203        off=i+j*d3_l->max_x+(k+1)*d3_l->max_x*d3_l->max_y;
204        carry=0;
205        if(!*(d3_l->status+off)&AMORPH) carry=(int)(my_info->dr_ac*(*(d3_l->extra+off)));
206        if(carry!=0)
207        {
208         *(d3_l->extra+off)-=carry;
209         *(d3_l->extra+offset)+=carry;
210        }
211       }
212      }  
213     }
214    } /* for z */
215   } /* for y */
216  } /* for x */
217
218  } /* if step modulo diff_rate == 0 */
219
220  return 1;
221 }
222
223 /* save to file --> send to server :)  --> T O D O <-- */
224 int send_data( ) {
225
226   return 1;
227 }
228
229 int save_to_file(char *sf,d3_lattice *d3_l,info *my_inf)
230 {
231  int sf_fd,c;
232
233  if((sf_fd=open(sf,O_WRONLY|O_CREAT))<0)
234  {
235   puts("cannot open save file");
236   return -1;
237  }
238  if(write(sf_fd,d3_l,sizeof(d3_lattice))<sizeof(d3_lattice))
239  {
240   puts("failed saving d3 lattice struct");
241   return -1;
242  }
243  if(write(sf_fd,my_inf,sizeof(info))<sizeof(info))
244  {
245   puts("failed saving info struct");
246   return-1;
247  }
248  c=d3_l->max_x*d3_l->max_y*d3_l->max_z;
249  if(write(sf_fd,d3_l->status,c*sizeof(unsigned char))<c*sizeof(unsigned char))
250  {
251   puts("failed saving status of d3 lattice sites");
252   return -1;
253  }
254  if(write(sf_fd,d3_l->extra,c*sizeof(int))<c*sizeof(int))
255  {
256   puts("failed saving sites concentration");
257   return -1;
258  }
259  close(sf_fd);
260
261  return 1;
262 }
263
264 /* load from file --> receive from server :) --> T O D O <-- */
265 int receive_data(  ) {
266
267   return 1;
268 }
269
270 int load_from_file(char *lf,d3_lattice *d3_l,info *my_inf)
271 {
272
273  int lf_fd,c,pos,end,data,data_len,strip;
274
275  if((lf_fd=open(lf,O_RDONLY))<0)
276  {
277   puts("cannot open load file");
278   return -1;
279  }
280  if(read(lf_fd,d3_l,sizeof(d3_lattice))<sizeof(d3_lattice))
281  {
282   puts("failed reading d3 lattice struct");
283   return -1;
284  }
285  pos=lseek(lf_fd,0,SEEK_CUR);
286  printf("psition: %d (%d)\n",pos,sizeof(d3_lattice));
287  data=d3_l->max_x*d3_l->max_y*d3_l->max_z;
288  data_len=data*(sizeof(int)+sizeof(unsigned char));
289  printf("there are %d volumes so we need %d of bytes\n",data,data_len);
290  end=lseek(lf_fd,0,SEEK_END);
291  c=end-pos-data_len;
292  printf("end: %d => length: %d => guessed info size: %d bytes\n",end,end-pos,c);
293  strip=sizeof(info)-c;
294  printf("as real programs info size is %d, we strip %d bytes\n",sizeof(info),strip);
295  lseek(lf_fd,pos,SEEK_SET);
296  c=sizeof(info);
297  if(strip>0) c-=strip;
298  if(c<0)
299  {
300   puts("info smaller then strip size");
301   return -1;
302  }
303  if(read(lf_fd,my_inf,c)<c)
304  {
305   puts("failed reading info struct");
306   return-1;
307  }
308  if(strip>0) memset(my_inf+c,0,strip);
309  if(strip<0) lseek(lf_fd,(-1*strip),SEEK_CUR);
310  c=d3_l->max_x*d3_l->max_y*d3_l->max_z;
311  if((d3_l->status=(unsigned char*)malloc(c*sizeof(unsigned char)))==NULL)
312  {
313   puts("cannot allocate status buffer");
314   return -1;
315  }
316  if((d3_l->extra=(int *)malloc(c*sizeof(int)))==NULL)
317  {
318   puts("cannot allocate concentration buffer");
319   return -1;
320  }
321  if(read(lf_fd,d3_l->status,c*sizeof(unsigned char))<c*sizeof(unsigned char))
322  {
323   puts("failed reading status of d3 lattice sites");
324   return -1;
325  }
326  if(read(lf_fd,d3_l->extra,c*sizeof(int))<c*sizeof(int))
327  {
328   puts("failed reading sites concentration");
329   return -1;
330  }
331  close(lf_fd);
332
333  return 1;
334 }
335
336 u32 get_reject_graph(info *my_info,d3_lattice *d3_l,char *file,u32 *graph) {
337  double a,b;
338  int i,j,k;
339  int fd;
340  char buf[32],*p;
341  unsigned char *flag;
342  u32 max;
343
344  max=0;
345  if((fd=open(file,O_RDONLY))<0)
346  {
347   puts("cannot open file to calculate rejection graph");
348   return -1;
349  }
350  if((flag=(unsigned char *)malloc(d3_l->max_z))==NULL)
351  {
352   puts("cannot malloc flag memory for rejection graph");
353   return -1;
354  }
355  memset(flag,0,d3_l->max_z);
356  memset(graph,0,d3_l->max_z*sizeof(u32));
357  /* get fixpoints */
358  k=1;
359  while(k)
360  {
361   for(i=0;i<32;i++)
362   {
363    k=read(fd,&buf[i],1);
364    if((buf[i]=='\n')||(k==0)) break;
365   }
366   if(k)
367   {
368    p=strtok(buf," ");
369    a=atof(p)/10; /* nm */
370    p=strtok(NULL," ");
371    b=atof(p);
372    if(a>d3_l->max_z*CELL_LENGTH) k=0;
373    else 
374    {
375     graph[(int)(a/CELL_LENGTH)]=(int)(URAND_MAX/100*b);
376     flag[(int)(a/CELL_LENGTH)]=1;
377    }
378   }
379  }
380  /* do (linear) interpolation here! */
381  i=0;
382  a=0;
383  while(i<d3_l->max_z)
384  {
385   /* graph[0] is 0! */
386   j=i;
387   i++;
388   while(flag[i]==0&&i<d3_l->max_z) i++;
389   for(k=j+1;k<i;k++) graph[k]=(int)((k-j)*((int)graph[i]-(int)graph[j])/(i-j))+graph[j];
390   if(graph[i]>max) max=graph[i];
391  }
392
393  free(flag);
394
395 #ifdef DEBUG_INTERPOL_PROFILE
396  printf("debug: %s (interpolated profile)\n",file);
397  for(i=0;i<d3_l->max_z;i++) printf("%d %d\n",i,graph[i]);
398 #endif
399
400  return max;
401 }
402
403
404
405 /*
406  * main program
407  */
408
409 int main(int argc,char **argv)
410 {
411
412   u32 x,y,z,x_c,y_c,z_c;
413   int i,j;
414   int resave;
415   int c_step;
416   char server_ip[16];
417   int port;
418   char p_file[MAX_CHARS];
419   char n_e_file[MAX_CHARS];
420   char r_file[MAX_CHARS];
421   d3_lattice d3_l;
422   info my_info;
423   u32 *c_profile;
424   u32 *n_e_loss;
425   u32 ne_max,ip_max;
426   u32 *nel_z;
427   unsigned char err_dc,dc;
428 #define DC_QUIT (1<<0)
429 #define DC_OK (1<<1)
430 #define DC_END (1<<2)
431   t_net net;
432
433   /* default values */
434   resave=RESAVE;
435   c_step=0;
436   strcpy(server_ip,"");
437   strcpy(p_file,IMP_PROFILE);
438   strcpy(n_e_file,NEL_PROFILE);
439   strcpy(r_file,"");
440   ne_max=0;
441   ip_max=0;
442   err_dc=0;
443   dc=0;
444   port=1025;
445
446   /* parse/check argv */
447   for(i=1;i<argc;i++) {
448     if(argv[i][0]=='-') {
449       switch(argv[i][1]) {
450         case 'h':
451           usage(argv[0]);
452           return -1;
453         case 'i':
454           strncpy(server_ip,argv[++i],16);
455           break;
456         case 'r':
457           strcpy(r_file,argv[++i]);
458           break;
459         case 'P':
460           strcpy(p_file,argv[++i]);
461           break;
462         case 'n':
463           strcpy(n_e_file,argv[++i]);
464           break;
465         case 'p':
466           port=atoi(argv[++i]);
467           break;
468         default:
469           usage(argv[0]);
470           return -1;
471       }
472     }
473   }
474   if(!strcmp(server_ip,"")) {
475     usage(argv[0]);
476     return -1;
477   }
478
479   /* connect to server and notify about free resource */
480   network_init(&net,1);
481   network_set_connection_info(&net,0,server_ip,port);
482   network_connect(&net,0);
483
484   /* wait for job */
485
486   /* get info (+data) */
487
488   /* care for signals */
489   err_dc=DC_QUIT;
490   signal(SIGTERM,send_data(&d3_l,&my_info,i,&err_dc));
491
492   /* rand init */
493   if(!strcmp(r_file,"")) rand_init(NULL);
494   else rand_init(r_file);
495
496   /* compute graphs for random number rejection method */
497   if((c_profile=(u32 *)malloc(d3_l.max_z*sizeof(unsigned int)))==NULL) {
498     puts("failed allocating memory for carbon profile graph");
499     return -1;
500   }
501   if((n_e_loss=(u32 *)malloc(d3_l.max_z*sizeof(unsigned int)))==NULL) {
502     puts("failed allocating memory for nuclear energy loss graph");
503     return -1;
504   }
505   ip_max=get_reject_graph(&my_info,&d3_l,p_file,c_profile);
506   ne_max=get_reject_graph(&my_info,&d3_l,n_e_file,n_e_loss);
507
508   /* array nel_z[] keeping nuclear energy loss values scaled to URAND_MAX */
509   nel_z=(u32 *)malloc(d3_l.max_z*sizeof(unsigned int));
510   if(nel_z==NULL) {
511     printf("failed allocating nel_z array mem\n");
512     return -1;
513   }
514   for(i=0;i<d3_l-max_z;i++) nel_z[i]=URAND_MAX*(1.0*n_e_loss[i]/ne_max);
515
516
517   /* this should be obsolete - z is high enough - we check now! */
518   if(c_profile[d3_l.max_z-1]!=0) {
519     printf("max_z (%d) too small - sputtering not possible\n",d3_l.max_z);
520     return -1;
521   }
522
523   /* sputtering really possible ?*/
524   if(n_e_loss[d3_l.max_z-1]!=0)
525     printf("warning: max_z (%d) too small, there may be amorphous volumes\n",
526            d3_l.max_z);
527
528   /*
529    * start of simulation
530    */
531
532   i=(c_step?c_step:0);
533   while(i<my_info.steps) {
534     for(j=0;j<my_info.cpi;j++) {
535       x_c=get_rand(d3_l.max_x);
536       y_c=get_rand(d3_l.max_y);
537       // z_c=get_rand_reject(d3_l.max_z,ne_max,n_e_loss);
538       z_c=get_rand(d3_l.max_z);
539       process_cell(&d3_l,x_c,y_c,z_c,&my_info,nel_z[z_c]);
540     }
541     distrib_c(&d3_l,&my_info,i,ip_max,c_profile);
542     if(i%resave==0 && i!=0) {
543       dc=DC_OK;
544       send_data(&d3_l,&my_info,i,&dc);
545     }
546     i++;
547     if(i%my_info.s_rate==0) sputter(&d3_l);
548   }
549
550   /* finished */
551   dc=DC_END;
552   send_data(&d3_l,&my_info,i,&dc);
553
554   /* shutdown/free/close everything now ... */
555
556   return 1;
557 }