missing ","
[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 <signal.h>
44
45 #include "nlsop.h"
46 #include "dfbapi.h"
47 #include "random.h"
48
49 #include "network.h"
50 #include "event.h"
51
52 #include "nlsop_general.h"
53
54 #define MAKE_AMORPH(N) *(N)|=AMORPH
55 #define MAKE_CRYST(N) *(N)&=~AMORPH
56
57 /* globals */
58
59 char p_file[MAX_CHARS];
60 char n_e_file[MAX_CHARS];
61 char r_file[MAX_CHARS];
62 t_net *gnet;
63 d3_lattice *gd3_l;
64 info *gmy_info;
65 int *gi;
66 unsigned char dc;
67
68 int get_data_and_calc(t_event *event,void *allineed);
69
70 int usage(char *prog)
71 {
72  puts("usage:");
73  printf("%s -i ip -p port -r/P/n random/profile/neloss file\n",prog);
74  return 1;
75 }
76
77 /*
78  * nlsop internal functions
79  */
80
81 int sputter(d3_lattice *d3_l)
82 {
83  int i,size;
84  int offh,offl;
85
86  size=d3_l->max_x*d3_l->max_y;
87  offl=0;
88  offh=size;
89
90  for(i=0;i<d3_l->max_z-1;i++)
91  {
92   memcpy(d3_l->status+offl,d3_l->status+offh,size);
93   memcpy(d3_l->extra+offl,d3_l->extra+offh,size*sizeof(int));
94   offl=offh;
95   offh+=size;
96  }
97  memset(d3_l->status+offl,0,size);
98  memset(d3_l->extra+offl,0,size);
99
100  return 1;
101 }
102
103 int process_cell(d3_lattice *d3_l,u32 x,u32 y,u32 z,info *my_info,u32 nel_z)
104 {
105  unsigned char *thiz;
106  int *conc;
107  int i,j;
108  int off;
109  double p,q;
110
111  thiz=d3_l->status+x+y*d3_l->max_x+z*d3_l->max_x*d3_l->max_y;
112  conc=d3_l->extra+x+y*d3_l->max_x+z*d3_l->max_x*d3_l->max_y;
113  p=my_info->b*nel_z;
114  for(i=-(my_info->range);i<=my_info->range;i++)
115  {
116   for(j=-(my_info->range);j<=my_info->range;j++)
117   {
118    if(!(i==0 && j==0))
119    {
120     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;
121     if(*(d3_l->status+off)&AMORPH) p+=my_info->s*(*(d3_l->extra+off))*URAND_MAX/(i*i+j*j);
122    } 
123   }
124  }
125  p+=*conc*my_info->c*URAND_MAX;
126  if(!(*thiz&AMORPH))
127  {
128   if(get_rand(URAND_MAX)<=p) MAKE_AMORPH(thiz);
129  } else
130  {
131   /* assume 1-p probability */
132   /* also look for neighbours ! */
133   q=(URAND_MAX-p)>0?URAND_MAX-p:0;
134   j=0;
135   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;
136   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;
137   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;
138   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;
139   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;
140   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;
141
142   p+=((q/6)*j);
143   if(get_rand(URAND_MAX)>p) MAKE_CRYST(thiz);
144  }
145  
146  return 1;
147 }
148
149 int distrib_c(d3_lattice *d3_l,info *my_info,int step,u32 rj_m,u32 *rj_g)
150 {
151  u32 x,y,z;
152  int i,j,k,c;
153  int offset,off;
154  int carry;
155
156  /* put one c ion somewhere in the lattice */
157  x=get_rand(d3_l->max_x);
158  y=get_rand(d3_l->max_y);
159  z=get_rand_reject(d3_l->max_z,rj_m,rj_g);
160  *(d3_l->extra+x+y*d3_l->max_x+z*d3_l->max_x*d3_l->max_y)+=1;
161  (my_info->cc)++;
162
163  if(step%my_info->diff_rate==0)
164  {
165
166  for(i=0;i<d3_l->max_x;i++)
167  {
168   for(j=0;j<d3_l->max_y;j++)
169   {
170    for(k=0;k<d3_l->max_z;k++)
171    {
172     offset=i+j*d3_l->max_x+k*d3_l->max_x*d3_l->max_y;
173     /* case amorph: amorph <- cryst diffusion */
174     if(*(d3_l->status+offset)&AMORPH)
175     {
176      for(c=-1;c<=1;c++)
177      {
178       if(c!=0)
179       {
180        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;
181        carry=0;
182        if(!(*(d3_l->status+off)&AMORPH)) carry=(int)(my_info->dr_ac*(*(d3_l->extra+off)));
183        if(carry!=0)
184        {
185         *(d3_l->extra+offset)+=carry;
186         *(d3_l->extra+off)-=carry;
187        }
188       }
189      }
190      for(c=-1;c<=1;c++)
191      {
192       if(c!=0)
193       {
194        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;
195        carry=0;
196        if(!(*(d3_l->status+off)&AMORPH)) carry=(int)(my_info->dr_ac*(*(d3_l->extra+off)));
197        if(carry!=0)
198        {
199         *(d3_l->extra+offset)+=carry; 
200         *(d3_l->extra+off)-=carry; 
201        }
202       }
203      }
204      /* diff in z direction */
205      if(k!=0)
206      {
207       off=i+j*d3_l->max_x+(k-1)*d3_l->max_x*d3_l->max_y;
208       carry=0;
209       if(!*(d3_l->status+off)&AMORPH) carry=(int)(my_info->dr_ac*(*(d3_l->extra+off)));
210       if(carry!=0)
211       {
212        *(d3_l->extra+off)-=carry;
213        *(d3_l->extra+offset)+=carry;
214       }
215      }
216      if(k!=d3_l->max_z-1)
217      {
218       off=i+j*d3_l->max_x+(k+1)*d3_l->max_x*d3_l->max_y;
219       carry=0;
220       if(!*(d3_l->status+off)&AMORPH) carry=(int)(my_info->dr_ac*(*(d3_l->extra+off)));
221       if(carry!=0)
222       {
223        *(d3_l->extra+off)-=carry;
224        *(d3_l->extra+offset)+=carry;
225       }
226      }
227     }
228    } /* for z */
229   } /* for y */
230  } /* for x */
231
232  } /* if step modulo diff_rate == 0 */
233
234  return 1;
235 }
236
237 u32 get_reject_graph(info *my_info,d3_lattice *d3_l,char *file,u32 *graph) {
238  double a,b;
239  int i,j,k;
240  int fd;
241  char buf[32],*p;
242  unsigned char *flag;
243  u32 max;
244
245  max=0;
246  if((fd=open(file,O_RDONLY))<0)
247  {
248   puts("cannot open file to calculate rejection graph");
249   return -1;
250  }
251  if((flag=(unsigned char *)malloc(d3_l->max_z))==NULL)
252  {
253   puts("cannot malloc flag memory for rejection graph");
254   return -1;
255  }
256  memset(flag,0,d3_l->max_z);
257  memset(graph,0,d3_l->max_z*sizeof(u32));
258  /* get fixpoints */
259  k=1;
260  while(k)
261  {
262   for(i=0;i<32;i++)
263   {
264    k=read(fd,&buf[i],1);
265    if((buf[i]=='\n')||(k==0)) break;
266   }
267   if(k)
268   {
269    p=strtok(buf," ");
270    a=atof(p)/10; /* nm */
271    p=strtok(NULL," ");
272    b=atof(p);
273    if(a>d3_l->max_z*CELL_LENGTH) k=0;
274    else 
275    {
276     graph[(int)(a/CELL_LENGTH)]=(int)(URAND_MAX/100*b);
277     flag[(int)(a/CELL_LENGTH)]=1;
278    }
279   }
280  }
281  /* do (linear) interpolation here! */
282  i=0;
283  a=0;
284  while(i<d3_l->max_z)
285  {
286   /* graph[0] is 0! */
287   j=i;
288   i++;
289   while(flag[i]==0&&i<d3_l->max_z) i++;
290   for(k=j+1;k<i;k++) graph[k]=(int)((k-j)*((int)graph[i]-(int)graph[j])/(i-j))+graph[j];
291   if(graph[i]>max) max=graph[i];
292  }
293
294  free(flag);
295
296 #ifdef DEBUG_INTERPOL_PROFILE
297  printf("debug: %s (interpolated profile)\n",file);
298  for(i=0;i<d3_l->max_z;i++) printf("%d %d\n",i,graph[i]);
299 #endif
300
301  return max;
302 }
303
304 void send_data(int signum) {
305
306   int c;
307
308   c=gd3_l->max_x*gd3_l->max_y*gd3_l->max_z;
309
310   network_send(gnet->connection[0].fd,&dc,1);
311   /*
312   network_send(gnet->connection[0].fd,(unsigned char *)gd3_l,
313                sizeof(d3_lattice));
314   network_send(gnet->connection[0].fd,(unsigned char *)gmy_info,sizeof(info));
315   */
316   network_send(gnet->connection[0].fd,gd3_l->status,c*sizeof(unsigned char));
317   network_send(gnet->connection[0].fd,(unsigned char *)gd3_l->extra,
318                c*sizeof(int));
319   network_send(gnet->connection[0].fd,(unsigned char *)gi,sizeof(int));
320
321 }
322
323
324 /*
325  * main program
326  */
327
328 int main(int argc,char **argv)
329 {
330
331   char server_ip[16];
332   int port;
333   t_net net;
334   t_event event;
335   unsigned char data[256];
336   int i;
337
338   gnet=&net;
339
340   /* default values */
341   strcpy(server_ip,"");
342   strcpy(p_file,IMP_PROFILE);
343   strcpy(n_e_file,NEL_PROFILE);
344   strcpy(r_file,"");
345   port=1025;
346
347   /* parse/check argv */
348   for(i=1;i<argc;i++) {
349     if(argv[i][0]=='-') {
350       switch(argv[i][1]) {
351         case 'h':
352           usage(argv[0]);
353           return -1;
354         case 'i':
355           strncpy(server_ip,argv[++i],16);
356           break;
357         case 'r':
358           strcpy(r_file,argv[++i]);
359           break;
360         case 'P':
361           strcpy(p_file,argv[++i]);
362           break;
363         case 'n':
364           strcpy(n_e_file,argv[++i]);
365           break;
366         case 'p':
367           port=atoi(argv[++i]);
368           break;
369         default:
370           usage(argv[0]);
371           return -1;
372       }
373     }
374   }
375   if(!strcmp(server_ip,"")) {
376     usage(argv[0]);
377     return -1;
378   }
379
380   /* event init */
381   event_init(&event,1);
382   event_set_timeout(&event,0,0);
383
384   /* connect to server */
385   network_init(&net,1);
386   network_set_connection_info(&net,0,server_ip,port);
387   if(network_connect(&net,0)==N_E_CONNECT) {
388     printf("unable to connect to server, aborting ...\n");
389     return -1;
390   }
391   network_select(&net,0);
392
393   /* tell server: i am a client, i may work for you */
394   data[0]=NLSOP_CLIENT;
395   network_send(net.connection[0].fd,data,1);
396
397   /* wait for job */
398   event_math(net.connection[0].fd,&event,READ,ADD);
399   printf("idle, waiting for jobs ...\n");
400   event_start(&event,NULL,get_data_and_calc,NULL);
401
402   return 1;
403 }
404
405 int get_data_and_calc(t_event *event,void *allineed) {
406
407   d3_lattice d3_l;
408   info my_info;
409   u32 *c_profile;
410   u32 *n_e_loss;
411   u32 ne_max,ip_max;
412   u32 *nel_z;
413   u32 x_c,y_c,z_c;
414   int i,j;
415   int c_step;
416   unsigned char data[256];
417
418   t_net *net;
419
420   c_step=0;
421   ne_max=0;
422   ip_max=0;
423
424   net=gnet;
425   gd3_l=&d3_l;
426   gmy_info=&my_info;
427   gi=&i;
428   dc=0;
429
430   printf("got a new job ...\n");
431   
432   /* get info (+data) */
433   network_receive(net->connection[0].fd,data,1);
434   if(data[0]==NLSOP_NJOB || data[0]==NLSOP_CJOB) {
435     network_receive(net->connection[0].fd,(unsigned char *)&d3_l,
436                     sizeof(d3_lattice));
437     network_receive(net->connection[0].fd,(unsigned char *)&my_info,
438                     sizeof(info));
439     c_step=0;
440     j=d3_l.max_x*d3_l.max_y*d3_l.max_z;
441     d3_l.status=(unsigned char *)malloc(j*sizeof(unsigned char));
442     if(d3_l.status==NULL) {
443       printf("status alloc failed\n");
444       return -1;
445     }
446     d3_l.extra=(int *)malloc(j*sizeof(int));
447     if(d3_l.extra==NULL) {
448       printf("extra malloc failed\n");
449       return -1;
450     }
451     if(data[0]==NLSOP_CJOB) {
452       network_receive(net->connection[0].fd,d3_l.status,
453                       j*sizeof(unsigned char));
454       network_receive(net->connection[0].fd,(unsigned char *)d3_l.extra,
455                       j*sizeof(int));
456       network_receive(net->connection[0].fd,(unsigned char *)&c_step,
457                       sizeof(int));
458     }
459   }
460
461   printf("starting simulation with following parameters:\n");
462   printf("b = %f | c = %f | s = %f\n",my_info.b,my_info.c,my_info.s);
463   printf("diff every %d steps | diff rate = %f\n",my_info.diff_rate,
464                                                   my_info.dr_ac);
465   printf("current step: %d | total steps: %d\n",my_info.steps,c_step);
466   printf("...\n");
467
468   /* care for signals */
469   dc=DC_QUIT;
470   signal(SIGTERM,send_data);
471
472   /* rand init */
473   if(!strcmp(r_file,"")) rand_init(NULL);
474   else rand_init(r_file);
475
476   /* compute graphs for random number rejection method */
477   if((c_profile=(u32 *)malloc(d3_l.max_z*sizeof(unsigned int)))==NULL) {
478     puts("failed allocating memory for carbon profile graph");
479     return -1;
480   }
481   if((n_e_loss=(u32 *)malloc(d3_l.max_z*sizeof(unsigned int)))==NULL) {
482     puts("failed allocating memory for nuclear energy loss graph");
483     return -1;
484   }
485   ip_max=get_reject_graph(&my_info,&d3_l,p_file,c_profile);
486   ne_max=get_reject_graph(&my_info,&d3_l,n_e_file,n_e_loss);
487
488   /* array nel_z[] keeping nuclear energy loss values scaled to URAND_MAX */
489   nel_z=(u32 *)malloc(d3_l.max_z*sizeof(unsigned int));
490   if(nel_z==NULL) {
491     printf("failed allocating nel_z array mem\n");
492     return -1;
493   }
494   for(i=0;i<d3_l.max_z;i++) nel_z[i]=URAND_MAX*(1.0*n_e_loss[i]/ne_max);
495
496   /* this should be obsolete - z is high enough - we check now! */
497   if(c_profile[d3_l.max_z-1]!=0) {
498     printf("max_z (%d) too small - sputtering not possible\n",d3_l.max_z);
499     return -1;
500   }
501
502   /* sputtering really possible ?*/
503   if(n_e_loss[d3_l.max_z-1]!=0)
504     printf("warning: max_z (%d) too small, there may be amorphous volumes\n",
505            d3_l.max_z);
506
507   /*
508    * start of simulation
509    */
510
511   i=(c_step?c_step:0);
512   while(i<my_info.steps) {
513     for(j=0;j<my_info.cpi;j++) {
514       x_c=get_rand(d3_l.max_x);
515       y_c=get_rand(d3_l.max_y);
516       // z_c=get_rand_reject(d3_l.max_z,ne_max,n_e_loss);
517       z_c=get_rand(d3_l.max_z);
518       process_cell(&d3_l,x_c,y_c,z_c,&my_info,nel_z[z_c]);
519     }
520     distrib_c(&d3_l,&my_info,i,ip_max,c_profile);
521     i++;
522     if(i%my_info.save_rate==0) {
523       dc=DC_OK;
524       send_data(0);
525       dc=DC_QUIT;
526     }
527     if(i%my_info.s_rate==0) sputter(&d3_l);
528   }
529
530   /* finished */
531   dc=DC_END;
532   send_data(0);
533   dc=DC_QUIT;
534
535   /* shutdown/free/close everything now ... */
536
537   return 1;
538 }