mknb_declarations.c

最著名最快的分子模拟软件

               
	if(mknb_func.coul==MKNB_COUL_TAB || mknb_func.coul==MKNB_COUL_GB || 
	   mknb_func.vdw==MKNB_VDW_TAB) {
		mknb_declare_real("r,rt,eps,eps2");
		mknb_declare_int("n0,nnn");
		mknb_declare_real("Y,F,Geps,Heps2,Fp,VV");
		if(mknb_func.do_force) {
			mknb_declare_real("FF");
			if(mknb_func.coul==MKNB_COUL_TAB || mknb_func.coul==MKNB_COUL_GB)
				mknb_declare_real("fijC");
			if(mknb_func.vdw==MKNB_VDW_TAB)
				mknb_declare_real("fijD,fijR");
		}
	}
	if(mknb_func.coul==MKNB_COUL_RF)
		mknb_declare_real("krsq");
	if(mknb_func.coul==MKNB_COUL_GB) {
		mknb_declare_real("isai,isaj,isaprod,gbscale");
		if(mknb_func.do_force) 
			mknb_declare_real("dvdasum,dvdatmp,dvdaj");
	}
    
	if(mknb_func.vdw==MKNB_VDW_BHAM) 
		mknb_declare_real("Vvdwexp,br");

	/* TIP4P water doesnt have any coulomb interaction
	 * on atom 1, so we skip it if we dont do LJ
	 */
	/* check for outer/i atom */
	firsti = ((mknb_func.vdw==MKNB_VDW_NO) &&
			  (mknb_func.water==MKNB_WATER_TIP4P_SINGLE || 
			   mknb_func.water==MKNB_WATER_TIP4P_PAIR)) ? 2 : 1;

	/* inner/j atom */
	firstj = ((mknb_func.vdw==MKNB_VDW_NO) && 
			  (mknb_func.water==MKNB_WATER_TIP4P_PAIR)) ? 2 : 1;
               
	/* i coordinates and forces */
	for(i=firsti;i<=mknb_func.ni;i++) {
		sprintf(buf,"ix%d,iy%d,iz%d",i,i,i);
		if(mknb_func.do_force) {
			sprintf(buf2,",fix%d,fiy%d,fiz%d",i,i,i);
			strcat(buf,buf2);
		}
		mknb_declare_real(buf);
	}
	/* j coordinates and forces */
	for(j=firstj;j<=mknb_func.nj;j++) {
		sprintf(buf,"jx%d,jy%d,jz%d",j,j,j);
		if(mknb_func.do_force && 
		   (mknb_func.water || mknb_options.prefetch_forces)) {
			if(!(mknb_func.water==MKNB_WATER_TIP4P_PAIR && j==1) || 
			   mknb_options.prefetch_forces) {
				sprintf(buf2,",fjx%d,fjy%d,fjz%d",j,j,j);
				strcat(buf,buf2);
			}
		}
		mknb_declare_real(buf);
	}
	/* i-j vectorial distance, rsq and rinv. */
	for(i=firsti;i<=mknb_func.ni;i++) {
		for(j=firstj;j<=mknb_func.nj;j++) {
			/* For TIP4p, site 1 never interacts with site 2,3,4 */
			if(mknb_func.water==MKNB_WATER_TIP4P_PAIR && 
			   ((i==1 && j>1) || (j==1 && i>1)))
				continue;
			sprintf(buf,"dx%d%d,dy%d%d,dz%d%d,rsq%d%d",
					i,j,i,j,i,j,i,j);
			if(mknb_func.coul || mknb_func.vdw!=MKNB_VDW_LJ) {
				if(!((mknb_func.water==MKNB_WATER_TIP4P_SINGLE || 
					  mknb_func.water==MKNB_WATER_TIP4P_PAIR) && 
					 i==1 && j==1 && mknb_func.vdw==MKNB_VDW_LJ)) {
					sprintf(buf2,",rinv%d%d",i,j);
					strcat(buf,buf2);
				}
			}
			mknb_declare_real(buf);
		}
	}
  
	/* The water charges and VdW parameters dont change,
	 * so we can determine them outside the mknb_func.
	 */
	if(mknb_func.water==MKNB_WATER_SPC_SINGLE)
		mknb_declare_real("qO,qH");
	else if(mknb_func.water==MKNB_WATER_TIP4P_SINGLE) 
		mknb_declare_real("qH,qM");
  
	if(mknb_func.water==MKNB_WATER_SPC_PAIR)
		mknb_declare_real("qO,qH,qqOO,qqOH,qqHH");
	if(mknb_func.water==MKNB_WATER_TIP4P_PAIR)
		mknb_declare_real("qH,qM,qqMM,qqMH,qqHH");

	if(mknb_func.vdw==MKNB_VDW_BHAM)
		mknb_declare_real("c6,cexp1,cexp2");
	else if(mknb_func.vdw!=MKNB_VDW_NO)
		mknb_declare_real("c6,c12");
  
	/* Variables needed for the inlined software inverse square root */
	if(mknb_options.software_invsqrt && 
	   (mknb_func.coul || mknb_func.vdw==MKNB_VDW_BHAM || 
		mknb_func.vdw==MKNB_VDW_TAB)) {
		mknb_declare_const_int("fractshift",12);
		mknb_declare_const_int("fractmask",8388607);
		mknb_declare_const_int("expshift",23);
		mknb_declare_const_int("expmask",2139095040);
		mknb_declare_const_int("explsb",8388608);
		mknb_declare_real4("lu");
		mknb_declare_int4("iexp,addr");
		/* To do bitwise manipulation of a FP number we need to move
		 * it back and forth between floating point and integer registers,
		 * without converting the actual data.
		 */
		if(mknb_fortran) {
			mknb_declare_int4("bval,result");
			mknb_declare_real4("fval");
			mknb_code("equivalence(bval,fval)");
			mknb_code("equivalence(result,lu)");
			mknb_declare_int4("invsqrtexptab,invsqrtfracttab");
			mknb_code("common /gmxinvsqrtdata/ invsqrtexptab(256),invsqrtfracttab(4096)");
		} else {
			mknb_declare_other("union { unsigned int bval; float fval; }",
						  "bitpattern,result");
		} 
	}
	fprintf(mknb_output,"\n");
}

void
mknb_initialize_data(void)
{   
	char buf[255];

	/* move arguments passed by reference to local scalars */
	if(!mknb_fortran) {
		mknb_assign("nri",        "*p_nri");
		mknb_assign("ntype",      "*p_ntype");
		mknb_assign("nthreads",   "*p_nthreads");
		mknb_assign("facel",      "*p_facel");
		mknb_assign("krf",        "*p_krf");
		mknb_assign("crf",        "*p_crf");
		mknb_assign("tabscale",   "*p_tabscale");
		if(mknb_func.coul==MKNB_COUL_GB)
    	{
			mknb_assign("gbtabscale", "*p_gbtabscale");
		}
	}	   

	/* assign the charge combinations for OO,OH and HH, 
     * or HH/HL/LL for TIP4P/TIP5P 
	 */
	/* we're always doing coulomb */
	if(mknb_func.water) {
		mknb_comment("Initialize water data");

		mknb_assign("ii", "%s%s", 
					mknb_array("iinr",(mknb_fortran) ? "1" : "0"), 
					(mknb_fortran) ? "+1" : "");
		switch(mknb_func.water) {
		case MKNB_WATER_SPC_SINGLE:
			mknb_assign("qO", "facel*%s", mknb_array("charge","ii"));
			mknb_assign("qH", "facel*%s", mknb_array("charge","ii+1"));
			break;
		case MKNB_WATER_TIP4P_SINGLE:
			mknb_assign("qH", "facel*%s", mknb_array("charge","ii+1"));
			mknb_assign("qM", "facel*%s", mknb_array("charge","ii+3"));
			break;
		case MKNB_WATER_SPC_PAIR:
			mknb_assign("qO", mknb_array("charge","ii"));
	 		mknb_assign("qH", mknb_array("charge","ii+1"));
			mknb_assign("qqOO","facel*qO*qO");
			mknb_assign("qqOH","facel*qO*qH");
			mknb_assign("qqHH","facel*qH*qH");
			break;
		case MKNB_WATER_TIP4P_PAIR:
			mknb_assign("qH", mknb_array("charge","ii+1"));
			mknb_assign("qM", mknb_array("charge","ii+3"));
			mknb_assign("qqMM","facel*qM*qM");
			mknb_assign("qqMH","facel*qM*qH");
			mknb_assign("qqHH","facel*qH*qH");
			break;
		default:
			printf("Error, unidentified water model (mknb_declarations.c)\n");
			exit(0);
		}

		if((mknb_func.water==MKNB_WATER_SPC_SINGLE || 
			mknb_func.water==MKNB_WATER_TIP4P_SINGLE) && mknb_func.vdw)
        {
			mknb_assign("nti","%d*ntype*%s",
                        mknb_func.nvdw_parameters,
                        mknb_array("type","ii"));
		}
        
		/* assign the nonbonded combination for the
		 * Oxygen-oxygen interactions 
		 */
		if((mknb_func.water==MKNB_WATER_SPC_PAIR || 
			mknb_func.water==MKNB_WATER_TIP4P_PAIR) && mknb_func.vdw) {
			sprintf(buf,"%d*(ntype+1)*%s%s",mknb_func.nvdw_parameters,
					mknb_array("type","ii"), (mknb_fortran) ? "+1" : "");
			mknb_assign("tj",buf);
			mknb_assign("c6",mknb_array("vdwparam","tj"));
			if(mknb_func.vdw==MKNB_VDW_BHAM) {
				mknb_assign("cexp1",mknb_array("vdwparam","tj+1"));
				mknb_assign("cexp2",mknb_array("vdwparam","tj+2"));
			} else
				mknb_assign("c12",mknb_array("vdwparam","tj+1"));
		}
		fprintf(mknb_output,"\n");
	}

	mknb_comment("Reset outer and inner iteration counters");

	if(mknb_options.threads)
	{
	    mknb_assign("nouter","0");
		mknb_assign("ninner","0");
	} 
	else
	{
		mknb_comment("Avoid compiler warning about unassigned variable");
		mknb_assign("nj1","0");
	}
}