From: hackbard <hackbard>
Date: Thu, 17 Aug 2006 15:41:43 +0000 (+0000)
Subject: more ..
X-Git-Url: https://hackdaworld.org/gitweb/?a=commitdiff_plain;h=e763f6e3072fbbbc5ecf7345b989d7d29929a38f;p=physik%2Fposic.git

more ..
---

diff --git a/moldyn.c b/moldyn.c
index 1b52a56..ad0f9a3 100644
--- a/moldyn.c
+++ b/moldyn.c
@@ -1005,6 +1005,9 @@ int tersoff(t_moldyn *moldyn) {
 				lambda=;
 				B=;
 				mu=;
+				chi=;
+				beta=;
+				betaN=;
 
 				if(d_ij<=R) {
 					f_c=1.0;
@@ -1032,6 +1035,10 @@ int tersoff(t_moldyn *moldyn) {
 			ck=link_cell_neighbour_index(moldyn,ki,kj,kk,
 			                             neighbourk);
 
+			/* go for zeta - 3 body stuff! */
+			zeta=0.0;
+			d_ij2=d_ij*d_ij;
+
 			/* cell of btom */
 			thisk=&(neighbourk[0]);
 			list_reset(thisk);
@@ -1044,12 +1051,45 @@ int tersoff(t_moldyn *moldyn) {
 				
 				/* 3 body stuff (1) */
 				
-				theta_ijk=;
-				sin_theta=;
-				cos_theta=;
-				hi_cos=;
-				hi_cos_square=;
+				v3_sub(&dist_ik,ktom,&(atom[i]));
+				d_ik=v3_norm(&dist_ik);
+				if(d_ik<=Sik) {
+
+				Rik=;
+				Sik=;
+				Aik=;
+				lambda_ik=;
+				Bik=;
+				mu_ik=;
+				omega_ik=;
+				c_i=;
+				d_i=;
+				h_i=;
+			
+
+					if(d_ik<=Rik) {
+						f_cik=1.0;
+						df_cik=0.0;
+					}
+					else {
+						sik_rik=Sik-Rik;
+						arg1ik=PI*(d_ik-Rik)/sik_rik;
+						f_cik=0.5+0.5*cos(arg1ik);
+						df_cik=-0.5*sin(arg1ik)* \
+						       (PI/(sik_rik*d_ik));
+						f_rik=Aik*exp(-lambda_ik*d_ik);
+						f_aik=-Bik*exp(-mu_ik*d_ik);
+					}
+			
+					v3_sub(&distance_jk,ktom,btom);
+					cos_theta=(d_ij2+d_ik*d_ik-d_jk*d_jk)/\
+					          (2*d_ij*d_ik);
+					theta=arccos(cos_theta);
 
+					
+				}
+				else
+					continue;
 
 				/* end 3 body stuff (1) */