X-Git-Url: https://www.hackdaworld.org/gitweb/?p=physik%2Fposic.git;a=blobdiff_plain;f=sic.c;h=375c6d54b38756aef74f12431075dd8fd6b446fc;hp=60df74b93b9e96318c11b3bc2056ef9cbf1b45d9;hb=HEAD;hpb=043b240b2312c1166440d11928e12fa3ad187086 diff --git a/sic.c b/sic.c index 60df74b..375c6d5 100644 --- a/sic.c +++ b/sic.c @@ -75,8 +75,8 @@ int insert_atoms(t_moldyn *moldyn) { // 001 dumbbell r.x=(-0.5+0.25)*ALBE_LC_SI; r.y=(-0.5+0.25)*ALBE_LC_SI; - r.z=(-0.5+0.25+0.125)*ALBE_LC_SI; - moldyn->atom[4372].r.z=(-0.5+0.25-0.125)*ALBE_LC_SI; + r.z=(-0.1)*ALBE_LC_SI; + moldyn->atom[4372].r.z=(-0.4)*ALBE_LC_SI; #endif #ifdef INS_USER // 001 dumbbell @@ -86,9 +86,18 @@ int insert_atoms(t_moldyn *moldyn) { #endif #ifdef INS_RAND // random +#ifdef INS_DYNAMIC_LEN + r.x=(rand_get_double(&(moldyn->random))-0.5)*\ + moldyn->dim.x; + r.y=(rand_get_double(&(moldyn->random))-0.5)*\ + moldyn->dim.y; + r.z=(rand_get_double(&(moldyn->random))-0.5)*\ + moldyn->dim.z; +#else r.x=(rand_get_double(&(moldyn->random))-0.5)*INS_LENX; r.y=(rand_get_double(&(moldyn->random))-0.5)*INS_LENY; r.z=(rand_get_double(&(moldyn->random))-0.5)*INS_LENZ; +#endif #endif // offset r.x+=INS_OFFSET; @@ -112,7 +121,7 @@ int insert_atoms(t_moldyn *moldyn) { dmin=d; } } - add_atom(moldyn,INS_TYPE,INS_MASS,INS_BRAND, + add_atom(moldyn,INS_TYPE,INS_BRAND, ATOM_ATTR_1BP|ATOM_ATTR_2BP|ATOM_ATTR_3BP|\ INS_ATTR, &r,&v); @@ -130,6 +139,7 @@ int sic_hook(void *moldyn,void *hook_params) { int steps; double tau; double dt; + double dp; hp=hook_params; md=moldyn; @@ -139,7 +149,8 @@ int sic_hook(void *moldyn,void *hook_params) { /* switch on t scaling */ if(md->schedule.count==0) - set_pt_scale(md,0,0,T_SCALE_BERENDSEN,T_SCALE_TAU); + set_pt_scale(md,P_SCALE_BERENDSEN,P_SCALE_TAU, + T_SCALE_BERENDSEN,T_SCALE_TAU); /* my lousy state machine ! */ @@ -163,24 +174,27 @@ int sic_hook(void *moldyn,void *hook_params) { insert: - /* immediately go on if no job is to be done */ - if(hp->insert_count==INS_RUNS) { - printf(" --- insert run return ---\n"); - hp->state=STATE_POSTRUN; - goto postrun; - } - /* assigne values */ steps=INS_RELAX; tau=INS_TAU; /* check temperature */ dt=md->t_avg-md->t_ref; + dp=md->p_avg-md->p_ref; if(dt<0) dt=-dt; - if(dt>INS_DELTA_TC) + if(dp<0) + dp=-dp; + if((dt>INS_DELTA_TC)|(dp>INS_DELTA_PC)) goto addsched; + /* immediately go on if no job is to be done */ + if(hp->insert_count==INS_RUNS) { + printf(" --- leaving insert state ---\n"); + hp->state=STATE_POSTRUN; + goto postrun; + } + /* else -> insert atoms */ hp->insert_count+=1; printf(" ### insert atoms (%d/%d) ###\n", @@ -190,23 +204,26 @@ insert: postrun: - /* immediately return if no job is to be done */ - if(hp->postrun_count==POST_RUNS) { - printf(" --- post run return ---\n"); - return 0; - } - /* assigne values */ steps=POST_RELAX; tau=POST_TAU; /* check temperature */ dt=md->t_avg-md->t_ref; + dp=md->p_avg-md->p_ref; if(dt<0) dt=-dt; - if(dt>POST_DELTA_TC) + if(dp<0) + dp=-dp; + if((dt>POST_DELTA_TC)|(dp>POST_DELTA_PC)) goto addsched; + /* immediately return if no job is to be done */ + if(hp->postrun_count==POST_RUNS) { + printf(" --- leaving post run state ---\n"); + return 0; + } + /* postrun action */ hp->postrun_count+=1; printf(" ### postrun (%d/%d) ###\n", @@ -232,10 +249,6 @@ int main(int argc,char **argv) { /* hook parameter structure */ t_hp hookparam; - /* potential parameters */ - t_tersoff_mult_params tp; - t_albe_mult_params ap; - /* testing location & velocity vector */ t_3dvec r,v; memset(&r,0,sizeof(t_3dvec)); @@ -249,32 +262,19 @@ int main(int argc,char **argv) { /* choose potential */ #ifdef ALBE - set_potential3b_j1(&md,albe_mult_3bp_j1); - set_potential3b_k1(&md,albe_mult_3bp_k1); - set_potential3b_j2(&md,albe_mult_3bp_j2); - set_potential3b_k2(&md,albe_mult_3bp_k2); + if(set_potential(&md,MOLDYN_POTENTIAL_AM)<0) + return -1; #else - set_potential1b(&md,tersoff_mult_1bp); - set_potential3b_j1(&md,tersoff_mult_3bp_j1); - set_potential3b_k1(&md,tersoff_mult_3bp_k1); - set_potential3b_j2(&md,tersoff_mult_3bp_j2); - set_potential3b_k2(&md,tersoff_mult_3bp_k2); -#endif - -#ifdef ALBE - set_potential_params(&md,&ap); -#else - set_potential_params(&md,&tp); + if(set_potential(&md,MOLDYN_POTENTIAL_TM)<0) + return -1; #endif /* cutoff radius & bondlen */ #ifdef ALBE set_cutoff(&md,ALBE_S_SI); - set_bondlen(&md,ALBE_S_SI,ALBE_S_C,ALBE_S_SIC); //set_cutoff(&md,ALBE_S_C); #else set_cutoff(&md,TM_S_SI); - set_bondlen(&md,TM_S_SI,TM_S_C,-1.0); //set_cutoff(&md,TM_S_C); #endif @@ -282,77 +282,17 @@ int main(int argc,char **argv) { * potential parameters */ +#ifndef ALBE /* * tersoff mult potential parameters for SiC */ - tp.S[0]=TM_S_SI; - tp.R[0]=TM_R_SI; - tp.A[0]=TM_A_SI; - tp.B[0]=TM_B_SI; - tp.lambda[0]=TM_LAMBDA_SI; - tp.mu[0]=TM_MU_SI; - tp.beta[0]=TM_BETA_SI; - tp.n[0]=TM_N_SI; - tp.c[0]=TM_C_SI; - tp.d[0]=TM_D_SI; - tp.h[0]=TM_H_SI; - - tp.S[1]=TM_S_C; - tp.R[1]=TM_R_C; - tp.A[1]=TM_A_C; - tp.B[1]=TM_B_C; - tp.lambda[1]=TM_LAMBDA_C; - tp.mu[1]=TM_MU_C; - tp.beta[1]=TM_BETA_C; - tp.n[1]=TM_N_C; - tp.c[1]=TM_C_C; - tp.d[1]=TM_D_C; - tp.h[1]=TM_H_C; - - tp.chi=TM_CHI_SIC; - - tersoff_mult_complete_params(&tp); - + tersoff_mult_set_params(&md,SI,C); +#else /* * albe mult potential parameters for SiC */ - ap.S[0]=ALBE_S_SI; - ap.R[0]=ALBE_R_SI; - ap.A[0]=ALBE_A_SI; - ap.B[0]=ALBE_B_SI; - ap.r0[0]=ALBE_R0_SI; - ap.lambda[0]=ALBE_LAMBDA_SI; - ap.mu[0]=ALBE_MU_SI; - ap.gamma[0]=ALBE_GAMMA_SI; - ap.c[0]=ALBE_C_SI; - ap.d[0]=ALBE_D_SI; - ap.h[0]=ALBE_H_SI; - - ap.S[1]=ALBE_S_C; - ap.R[1]=ALBE_R_C; - ap.A[1]=ALBE_A_C; - ap.B[1]=ALBE_B_C; - ap.r0[1]=ALBE_R0_C; - ap.lambda[1]=ALBE_LAMBDA_C; - ap.mu[1]=ALBE_MU_C; - ap.gamma[1]=ALBE_GAMMA_C; - ap.c[1]=ALBE_C_C; - ap.d[1]=ALBE_D_C; - ap.h[1]=ALBE_H_C; - - ap.Smixed=ALBE_S_SIC; - ap.Rmixed=ALBE_R_SIC; - ap.Amixed=ALBE_A_SIC; - ap.Bmixed=ALBE_B_SIC; - ap.r0_mixed=ALBE_R0_SIC; - ap.lambda_m=ALBE_LAMBDA_SIC; - ap.mu_m=ALBE_MU_SIC; - ap.gamma_m=ALBE_GAMMA_SIC; - ap.c_mixed=ALBE_C_SIC; - ap.d_mixed=ALBE_D_SIC; - ap.h_mixed=ALBE_H_SIC; - - albe_mult_complete_params(&ap); + albe_mult_set_params(&md,SI,C); +#endif /* set (initial) dimensions of simulation volume */ #ifdef ALBE @@ -385,25 +325,25 @@ int main(int argc,char **argv) { // diamond #ifdef ALBE #ifdef INIT_SI - create_lattice(&md,DIAMOND,ALBE_LC_SI,SI,M_SI, + create_lattice(&md,DIAMOND,ALBE_LC_SI,SI, ATOM_ATTR_1BP|ATOM_ATTR_2BP|ATOM_ATTR_3BP|ATOM_ATTR_HB, - 0,LCNTX,LCNTY,LCNTZ,NULL); + 0,LCNTX,LCNTY,LCNTZ,NULL,0,NULL); #endif #ifdef INIT_C - create_lattice(&md,DIAMOND,ALBE_LC_C,C,M_C, + create_lattice(&md,DIAMOND,ALBE_LC_C,C, ATOM_ATTR_1BP|ATOM_ATTR_2BP|ATOM_ATTR_3BP|ATOM_ATTR_HB, - 1,LCNTX,LCNTY,LCNTZ,NULL); + 1,LCNTX,LCNTY,LCNTZ,NULL,0,NULL); #endif #else #ifdef INIT_SI - create_lattice(&md,DIAMOND,LC_SI,SI,M_SI, + create_lattice(&md,DIAMOND,LC_SI,SI, ATOM_ATTR_1BP|ATOM_ATTR_2BP|ATOM_ATTR_3BP|ATOM_ATTR_HB, - 0,LCNTX,LCNTY,LCNTZ,NULL); + 0,LCNTX,LCNTY,LCNTZ,NULL,0,NULL); #endif #ifdef INIT_C - create_lattice(&md,DIAMOND,LC_C,SI,M_SI, + create_lattice(&md,DIAMOND,LC_C,SI, ATOM_ATTR_1BP|ATOM_ATTR_2BP|ATOM_ATTR_3BP|ATOM_ATTR_HB, - 1,LCNTX,LCNTY,LCNTZ,NULL); + 1,LCNTX,LCNTY,LCNTZ,NULL,0,NULL); #endif #endif @@ -411,22 +351,22 @@ int main(int argc,char **argv) { #ifdef INIT_3CSIC #ifdef ALBE r.x=0.5*0.25*ALBE_LC_SIC; r.y=r.x; r.z=r.x; - create_lattice(&md,FCC,ALBE_LC_SIC,SI,M_SI, + create_lattice(&md,FCC,ALBE_LC_SIC,SI, ATOM_ATTR_1BP|ATOM_ATTR_2BP|ATOM_ATTR_3BP|ATOM_ATTR_HB, - 0,LCNTX,LCNTY,LCNTZ,&r); + 0,LCNTX,LCNTY,LCNTZ,&r,0,NULL); r.x+=0.25*ALBE_LC_SIC; r.y=r.x; r.z=r.x; - create_lattice(&md,FCC,ALBE_LC_SIC,C,M_C, - ATOM_ATTR_1BP|ATOM_ATTR_2BP|ATOM_ATTR_3BP|ATOM_ATTR_HB, - 1,LCNTX,LCNTY,LCNTZ,&r); + create_lattice(&md,FCC,ALBE_LC_SIC,C, + ATOM_ATTR_1BP|ATOM_ATTR_2BP|ATOM_ATTR_3BP|ATOM_ATTR_HB|ATOM_ATTR_VB, + 1,LCNTX,LCNTY,LCNTZ,&r,0,NULL); #else r.x=0.5*0.25*TM_LC_SIC; r.y=r.x; r.z=r.x; - create_lattice(&md,FCC,TM_LC_SIC,SI,M_SI, + create_lattice(&md,FCC,TM_LC_SIC,SI, ATOM_ATTR_1BP|ATOM_ATTR_2BP|ATOM_ATTR_3BP|ATOM_ATTR_HB, - 0,LCNTX,LCNTY,LCNTZ,&r); + 0,LCNTX,LCNTY,LCNTZ,&r,0,NULL); r.x+=0.25*TM_LC_SIC; r.y=r.x; r.z=r.x; - create_lattice(&md,FCC,TM_LC_SIC,C,M_C, + create_lattice(&md,FCC,TM_LC_SIC,C, ATOM_ATTR_1BP|ATOM_ATTR_2BP|ATOM_ATTR_3BP|ATOM_ATTR_HB, - 1,LCNTX,LCNTY,LCNTZ,&r); + 1,LCNTX,LCNTY,LCNTZ,&r,0,NULL); #endif #endif @@ -467,14 +407,14 @@ int main(int argc,char **argv) { /* set temperature & pressure */ set_temperature(&md,atof(argv[2])+273.0); - set_pressure(&md,BAR); + set_pressure(&md,0.0); /* set amount of steps to skip before average calc */ set_avg_skip(&md,AVG_SKIP); /* set p/t scaling */ //set_pt_scale(&md,0,0,T_SCALE_BERENDSEN,100.0); - //set_pt_scale(&md,P_SCALE_BERENDSEN,0.001, + //set_pt_scale(&md,P_SCALE_BERENDSEN,0.01/(100*GPA), // T_SCALE_BERENDSEN,100.0); //set_pt_scale(&md,0,0,T_SCALE_DIRECT,1.0); //set_pt_scale(&md,P_SCALE_BERENDSEN,0.001,0,0); @@ -499,7 +439,8 @@ int main(int argc,char **argv) { moldyn_set_log(&md,LOG_TOTAL_ENERGY,LOG_E); moldyn_set_log(&md,LOG_TEMPERATURE,LOG_T); moldyn_set_log(&md,LOG_PRESSURE,LOG_P); - moldyn_set_log(&md,VISUAL_STEP,LOG_V); + moldyn_set_log(&md,LOG_VOLUME,LOG_V); + moldyn_set_log(&md,VISUAL_STEP,LOG_A); moldyn_set_log(&md,SAVE_STEP,LOG_S); moldyn_set_log(&md,CREATE_REPORT,0);