added harmonic potntial + bugfixes + boundings
[physik/posic.git] / moldyn.h
index d94e193..a2633ba 100644 (file)
--- a/moldyn.h
+++ b/moldyn.h
@@ -20,6 +20,7 @@ typedef struct s_atom {
        t_3dvec f;      /* forces */
        int element;    /* number of element in pse */
        double mass;    /* atom mass */
+       //t_list vicinity       /* verlet neighbour list */
 } t_atom;
 
 typedef struct s_moldyn {
@@ -28,15 +29,22 @@ typedef struct s_moldyn {
        double (*potential)(struct s_moldyn *moldyn);
        void *pot_params;
        int (*force)(struct s_moldyn *moldyn);
+       double cutoff;
        double cutoff_square;
        t_3dvec dim;
        int (*integrate)(struct s_moldyn *moldyn);
        int time_steps;
        double tau;
        void *visual;
+       int write;
        unsigned char status;
 } t_moldyn;
 
+typedef struct s_ho_params {
+       double spring_constant;
+       double equilibrium_distance;
+} t_ho_params;
+
 typedef struct s_lj_params {
        double sigma6;
        double sigma12;
@@ -54,18 +62,20 @@ typedef struct s_lj_params {
 
 /* phsical values */
 
-//#define K_BOLTZMANN          1.3807E-23
-#define K_BOLTZMANN            1.0
+#define K_BOLTZMANN            1.3807e-27                      /* Nm/K */
+#define AMU                    1.660540e-27                    /* kg */
 
 #define FCC                    0x01
 #define DIAMOND                        0x02
 
 #define C                      0x06
-#define M_C                    6.0
+#define M_C                    (12.011*AMU)
 
-#define Si                     0x0e
-#define M_SI                   14.0
-#define LC_SI                  5.43105
+#define SI                     0x0e
+#define LC_SI                  0.543105e-9                             /* m */
+#define M_SI                   (28.085*AMU)                            /* kg */
+#define LJ_SIGMA_SI            ((0.25*sqrt(3.0)*LC_SI)/1.122462)       /* m */
+#define LJ_EPSILON_SI          (2.1678*1.60e-19)                       /* Nm */
 
 /* function prototypes */
 
@@ -79,9 +89,13 @@ double get_e_pot(t_moldyn *moldyn);
 double get_total_energy(t_moldyn *moldyn);
 t_3dvec get_total_p(t_atom *atom,int count);
 
+double estimate_time_step(t_moldyn *moldyn,double nn_dist,double t);
+
 int moldyn_integrate(t_moldyn *moldyn);
 int velocity_verlet(t_moldyn *moldyn);
 
+double potential_harmonic_oscillator(t_moldyn *moldyn);
+int force_harmonic_oscillator(t_moldyn *moldyn);
 double potential_lennard_jones(t_moldyn *moldyn);
 int force_lennard_jones(t_moldyn *moldyn);