1 /*
2 *
3 * CDDL HEADER START
4 *
5 * The contents of this file are subject to the terms of the Common Development
6 * and Distribution License Version 1.0 (the "License").
7 *
8 * You can obtain a copy of the license at
9 * http://www.opensource.org/licenses/CDDL-1.0.  See the License for the
10 * specific language governing permissions and limitations under the License.
11 *
12 * When distributing Covered Code, include this CDDL HEADER in each file and
13 * include the License file in a prominent location with the name LICENSE.CDDL.
14 * If applicable, add the following below this CDDL HEADER, with the fields
15 * enclosed by brackets "[]" replaced with your own identifying information:
16 *
17 * Portions Copyright (c) [yyyy] [name of copyright owner]. All rights reserved.
18 *
19 * CDDL HEADER END
20 *
21 
22 *
23 * Copyright (c) 2012, Regents of the University of Minnesota.  All rights reserved.
24 *
25 * Contributors:
26 *    Ryan S. Elliott
27 *    Ellad B. Tadmor
28 *    Valeriu Smirichinski
29 *
30 */
31 
32 /*******************************************************************************
33 *
34 *  model_driver_P_<FILL model driver name>
35 *
36 *  <FILL model drive name> pair potential KIM Model Driver
37 *
38 *  Reference: <FILL>
39 *
40 *  Language: C
41 *
42 *  Release: This file is part of the openkim-api-v1.1.1 package.
43 *
44 *******************************************************************************/
45 
46 
47 #include <stdlib.h>
48 #include <stdio.h>
49 #include <string.h>
50 #include <math.h>
51 #include "KIM_API_C.h"
52 #include "KIM_API_status.h"
53 
54 /******************************************************************************
55 * Below are the definitions and values of all Model parameters
56 *******************************************************************************/
57 #define DIM 3       /* dimensionality of space */
58 #define SPECCODE 1  /* internal species code */
59 
60 
61 /* Define prototypes for Model Driver init */
62 /* must be all lowercase to be compatible with the KIM API (to support Fortran Tests) */
63 /**/
64 int model_driver_p_<FILL (lowercase) model driver name>_init_(void* km, char* paramfile_names, int* nmstrlen, int* numparamfiles);
65 
66 /* Define prototypes for Model (Driver) reinit, compute, and destroy */
67 /* defined as static to avoid namespace clashes with other Models    */
68 /**/
69 static int reinit(void* km);
70 static int destroy(void* km);
71 static int compute(void* km);
72 /**/
73 static void calc_phi(double* <FILL parameter 1>,
74                      double* <FILL parameter 2>,
75                      /* FILL as many parameters as needed */
76                      double* cutoff, double r, double* phi);
77 static void calc_phi_dphi(double* <FILL parameter 1>,
78                           double* <FILL parameter 2>,
79                           /* FILL as many parameters as needed */
80                           double* cutoff, double r, double* phi, double* dphi);
81 
82 static void calc_phi_d2phi(double* <FILL parameter 1>,
83                            double* <FILL parameter 2>,
84                            /* FILL as many parameters as needed */
85                            double* cutoff, double r, double* phi, double* dphi, double* d2phi);
86 
87 /* Define model_buffer structure */
88 struct model_buffer {
89    int NBC;
90    int HalfOrFull;
91    int IterOrLoca;
92    int energy_ind;
93    int forces_ind;
94    int particleEnergy_ind;
95    int process_dEdr_ind;
96    int process_d2Edr2_ind;
97    int model_index_shift;
98    int numberOfParticles_ind;
99    int particleTypes_ind;
100    int coordinates_ind;
101    int numberContributingParticles_ind;
102    int boxSideLengths_ind;
103    int get_neigh_ind;
104    int cutoff_ind;
105 
106 
107    double Pcutoff;
108    double cutsq;
109    double <FILL parameter 1>;
110    double <FILL parameter 2>;
111    /* FILL as many parameters as needed */
112 };
113 
114 
115 /* Calculate pair potential phi(r) */
calc_phi(double * <FILL parameter1>,double * <FILL parameter2>,double * cutoff,double r,double * phi)116 static void calc_phi(double* <FILL parameter 1>,
117                      double* <FILL parameter 2>,
118                      /* FILL as many parameters as needed */
119                      double* cutoff, double r, double* phi)
120 {
121    /* local variables */
122    /* FILL: place any local variable definitions here */
123 
124    if (r > *cutoff)
125    {
126       /* Argument exceeds cutoff radius */
127       *phi = 0.0;
128    }
129    else
130    {
131       *phi   = <FILL functional form of phi(r)>;
132    }
133 
134    return;
135 }
136 
137 /* Calculate pair potential phi(r) and its derivative dphi(r) */
calc_phi_dphi(double * <FILL parameter1>,double * <FILL parameter2>,double * cutoff,double r,double * phi,double * dphi)138 static void calc_phi_dphi(double* <FILL parameter 1>,
139                           double* <FILL parameter 2>,
140                           /* FILL as many parameters as needed */
141                           double* cutoff, double r, double* phi, double* dphi)
142 {
143    /* local variables */
144    /* FILL: place any local variable definitions here */
145 
146    if (r > *cutoff)
147    {
148       /* Argument exceeds cutoff radius */
149       *phi  = 0.0;
150       *dphi = 0.0;
151    }
152    else
153    {
154       *phi  = <FILL functional form of phi(r)>;
155       *dphi = <FILL functional form of dphi(r)>;
156    }
157 
158    return;
159 }
160 
161 /* Calculate pair potential phi(r) and its 1st & 2nd derivatives dphi(r), d2phi(r) */
calc_phi_d2phi(double * <FILL parameter1>,double * <FILL parameter2>,double * cutoff,double r,double * phi,double * dphi,double * d2phi)162 static void calc_phi_d2phi(double* <FILL parameter 1>,
163                            double* <FILL parameter 2>,
164                            /* FILL as many parameters as needed */
165                            double* cutoff, double r, double* phi, double* dphi, double* d2phi)
166 {
167    /* local variables */
168    /* FILL: place any local variable definitions here */
169 
170    if (r > *cutoff)
171    {
172       /* Argument exceeds cutoff radius */
173       *phi   = 0.0;
174       *dphi  = 0.0;
175       *d2phi = 0.0;
176    }
177    else
178    {
179       *phi   = <FILL functional form of phi(r)>;
180       *dphi  = <FILL functional form of dphi(r)>;
181       *d2phi = <FILL functional form of d2phi(r)>;
182    }
183 
184    return;
185 }
186 
187 /* compute function */
compute(void * km)188 static int compute(void* km)
189 {
190    /* local variables */
191    intptr_t* pkim = *((intptr_t**) km);
192    double R;
193    double R_pairs[2];
194    double *pR_pairs = &(R_pairs[0]);
195    double Rsqij;
196    double phi;
197    double dphi;
198    double d2phi;
199    double dEidr;
200    double d2Eidr;
201    double Rij[DIM];
202    double *pRij = &(Rij[0]);
203    double Rij_pairs[2][3];
204    double *pRij_pairs = &(Rij_pairs[0][0]);
205    int ier;
206    int i;
207    int i_pairs[2];
208    int *pi_pairs = &(i_pairs[0]);
209    int j;
210    int j_pairs[2];
211    int *pj_pairs = &(j_pairs[0]);
212    int jj;
213    int k;
214    int currentAtom;
215    int* neighListOfCurrentAtom;
216    struct model_buffer* buffer;
217    int comp_energy;
218    int comp_force;
219    int comp_particleEnergy;
220    int comp_process_dEdr;
221    int comp_process_d2Edr2;
222    int NBC;
223    int HalfOrFull;
224    int IterOrLoca;
225    int model_index_shift;
226    int zero = 0;
227    int one = 1;
228    int request;
229 
230    int* nAtoms;
231    int* particleTypes;
232    double* cutoff;
233    double* cutsq;
234    double* <FILL parameter 1>;
235    double* <FILL parameter 2>;
236    /* FILL as many parameters as needed */
237    double* Rij_list;
238    double* coords;
239    double* energy;
240    double* force;
241    double* particleEnergy;
242    double* boxSideLengths;
243    int* numContrib;
244    int numberContrib;
245    int numOfAtomNeigh;
246    int (*get_neigh)(void *,int *,int *,int *, int *, int **, double **);
247 
248 
249    /* get buffer from KIM object */
250    buffer = (struct model_buffer*) KIM_API_get_model_buffer(pkim, &ier);
251    if (KIM_STATUS_OK > ier)
252    {
253       KIM_API_report_error(__LINE__, __FILE__, "KIM_API_get_model_buffer", ier);
254       return ier;
255    }
256 
257    /* unpack info from the buffer */
258    NBC = buffer->NBC;
259    HalfOrFull = buffer->HalfOrFull;
260    IterOrLoca = buffer->IterOrLoca;
261    model_index_shift = buffer->model_index_shift;
262    cutsq = &(buffer->cutsq);
263    <FILL parameter 1> = &(buffer-><FILL parameter 1>);
264    <FILL parameter 2> = &(buffer-><FILL parameter 2>);
265    /* also FILL additional parameters here if there are any ... */
266 
267    /* check to see if we have been asked to compute the forces, particleEnergy, and d1Edr */
268    KIM_API_getm_compute_by_index(pkim, &ier, 5*3,
269                                  buffer->energy_ind,         &comp_energy,         1,
270                                  buffer->forces_ind,         &comp_force,          1,
271                                  buffer->particleEnergy_ind, &comp_particleEnergy, 1,
272                                  buffer->process_dEdr_ind,   &comp_process_dEdr,   1,
273                                  buffer->process_d2Edr2_ind, &comp_process_d2Edr2, 1);
274    if (KIM_STATUS_OK > ier)
275    {
276       KIM_API_report_error(__LINE__, __FILE__, "KIM_API_getm_compute_by_index", ier);
277       return ier;
278    }
279 
280    KIM_API_getm_data_by_index(pkim, &ier, 10*3,
281                               buffer->cutoff_ind,                      &cutoff,         1,
282                               buffer->numberOfParticles_ind,           &nAtoms,         1,
283                               buffer->particleTypes_ind,               &particleTypes,  1,
284                               buffer->coordinates_ind,                 &coords,         1,
285                               buffer->numberContributingParticles_ind, &numContrib,     (HalfOrFull==1),
286                               buffer->get_neigh_ind,                   &get_neigh,      (NBC!=0),
287                               buffer->boxSideLengths_ind,              &boxSideLengths, (NBC==1),
288                               buffer->energy_ind,                      &energy,         comp_energy,
289                               buffer->forces_ind,                      &force,          comp_force,
290                               buffer->particleEnergy_ind,              &particleEnergy, comp_particleEnergy);
291    if (KIM_STATUS_OK > ier)
292    {
293       KIM_API_report_error(__LINE__, __FILE__, "KIM_API_getm_data_by_index", ier);
294       return ier;
295    }
296 
297    if (HalfOrFull == 1)
298    {
299       if (0 != NBC) /* non-CLUSTER cases */
300       {
301          numberContrib = *numContrib;
302       }
303       else /* CLUSTER cases */
304       {
305          numberContrib = *nAtoms;
306       }
307    }
308 
309    /* Check to be sure that the atom types are correct */
310    /**/
311    ier = KIM_STATUS_FAIL; /* assume an error */
312    for (i = 0; i < *nAtoms; ++i)
313    {
314       if ( SPECCODE != particleTypes[i])
315       {
316          KIM_API_report_error(__LINE__, __FILE__, "Unexpected species type detected", ier);
317          return ier;
318       }
319    }
320    ier = KIM_STATUS_OK; /* everything is ok */
321 
322    /* initialize potential energies, forces, and virial term */
323    if (comp_particleEnergy)
324    {
325       for (i = 0; i < *nAtoms; ++i)
326       {
327          particleEnergy[i] = 0.0;
328       }
329    }
330    if (comp_energy)
331    {
332       *energy = 0.0;
333    }
334 
335    if (comp_force)
336    {
337       for (i = 0; i < *nAtoms; ++i)
338       {
339          for (k = 0; k < DIM; ++k)
340          {
341             force[i*DIM + k] = 0.0;
342          }
343       }
344    }
345 
346    /* Initialize neighbor handling for CLUSTER NBC */
347    if (0 == NBC) /* CLUSTER */
348    {
349       neighListOfCurrentAtom = (int *) malloc((*nAtoms)*sizeof(int));
350    }
351 
352    /* Initialize neighbor handling for Iterator mode */
353 
354    if (1 == IterOrLoca)
355    {
356       ier = (*get_neigh)(&pkim, &zero, &zero, &currentAtom, &numOfAtomNeigh,
357                           &neighListOfCurrentAtom, &Rij_list);
358       /* check for successful initialization */
359       if (KIM_STATUS_NEIGH_ITER_INIT_OK != ier)
360       {
361          KIM_API_report_error(__LINE__, __FILE__, "KIM_API_get_neigh", ier);
362          ier = KIM_STATUS_FAIL;
363          return ier;
364       }
365    }
366 
367    /* Compute energy and forces */
368 
369    /* loop over particles and compute enregy and forces */
370    i = -1;
371    while( 1 )
372    {
373 
374       /* Set up neighbor list for next atom for all NBC methods */
375       if (1 == IterOrLoca) /* ITERATOR mode */
376       {
377          ier = (*get_neigh)(&pkim, &zero, &one, &currentAtom, &numOfAtomNeigh,
378                              &neighListOfCurrentAtom, &Rij_list);
379          if (KIM_STATUS_NEIGH_ITER_PAST_END == ier) /* the end of the list, terminate loop */
380          {
381             break;
382          }
383          if (KIM_STATUS_OK > ier) /* some sort of problem, exit */
384          {
385             KIM_API_report_error(__LINE__, __FILE__, "KIM_API_get_neigh", ier);
386             return ier;
387          }
388 
389          i = currentAtom + model_index_shift;
390       }
391       else
392       {
393          i++;
394          if (*nAtoms <= i) /* incremented past end of list, terminate loop */
395          {
396             break;
397          }
398 
399          if (0 == NBC) /* CLUSTER NBC method */
400          {
401             numOfAtomNeigh = *nAtoms - (i + 1);
402             for (k = 0; k < numOfAtomNeigh; ++k)
403             {
404                neighListOfCurrentAtom[k] = i + k + 1 - model_index_shift;
405             }
406             ier = KIM_STATUS_OK;
407          }
408          else
409          {
410             request = i - model_index_shift;
411             ier = (*get_neigh)(&pkim, &one, &request,
412                                 &currentAtom, &numOfAtomNeigh,
413                                 &neighListOfCurrentAtom, &Rij_list);
414             if (KIM_STATUS_OK != ier) /* some sort of problem, exit */
415             {
416                KIM_API_report_error(__LINE__, __FILE__, "get_neigh", ier);
417                ier = KIM_STATUS_FAIL;
418                return ier;
419             }
420          }
421       }
422 
423       /* loop over the neighbors of atom i */
424       for (jj = 0; jj < numOfAtomNeigh; ++ jj)
425       {
426 
427          j = neighListOfCurrentAtom[jj] + model_index_shift; /* get neighbor ID */
428 
429          /* compute relative position vector and squared distance */
430          Rsqij = 0.0;
431          for (k = 0; k < DIM; ++k)
432          {
433             if (3 != NBC) /* all methods except NEIGH_RVEC_F */
434             {
435                Rij[k] = coords[j*DIM + k] - coords[i*DIM + k];
436             }
437             else          /* NEIGH_RVEC_F method */
438             {
439                Rij[k] = Rij_list[jj*DIM + k];
440             }
441 
442             /* apply periodic boundary conditions if required */
443             if (1 == NBC)
444             {
445                if (abs(Rij[k]) > 0.5*boxSideLengths[k])
446                {
447                   Rij[k] -= (Rij[k]/fabs(Rij[k]))*boxSideLengths[k];
448                }
449             }
450 
451             /* compute squared distance */
452             Rsqij += Rij[k]*Rij[k];
453          }
454 
455          /* compute energy and force */
456          if (Rsqij < *cutsq) /* particles are interacting ? */
457          {
458             R = sqrt(Rsqij);
459             if (comp_process_d2Edr2)
460             {
461                /* compute pair potential and its derivatives */
462                calc_phi_d2phi(<FILL parameter 1>,
463                               <FILL parameter 2>,
464                               /* FILL as many parameters as needed */
465                               cutoff, R, &phi, &dphi, &d2phi);
466 
467                /* compute dEidr */
468                if ((1 == HalfOrFull) && (j < numberContrib))
469                {
470                   /* Half mode -- double contribution */
471                   dEidr = dphi;
472                   d2Eidr = d2phi;
473                }
474                else
475                {
476                   /* Full mode -- regular contribution */
477                   dEidr = 0.5*dphi;
478                   d2Eidr = 0.5*d2phi;
479                }
480             }
481             else if (comp_force || comp_process_dEdr)
482             {
483                /* compute pair potential and its derivative */
484                calc_phi_dphi(<FILL parameter 1>,
485                              <FILL parameter 2>,
486                              /* FILL as many parameters as needed */
487                              cutoff, R, &phi, &dphi);
488 
489                /* compute dEidr */
490                if ((1 == HalfOrFull) && (j < numberContrib))
491                {
492                   /* Half mode -- double contribution */
493                   dEidr = dphi;
494                }
495                else
496                {
497                   /* Full mode -- regular contribution */
498                   dEidr = 0.5*dphi;
499                }
500             }
501             else
502             {
503                /* compute just pair potential */
504                calc_phi(<FILL parameter 1>,
505                         <FILL parameter 2>,
506                         /* FILL as many parameters as needed */
507                         cutoff, R, &phi);
508             }
509 
510             /* contribution to energy */
511             if (comp_particleEnergy)
512             {
513                particleEnergy[i] += 0.5*phi;
514                /* if half list add energy for the other atom in the pair */
515                if ((1 == HalfOrFull) && (j < numberContrib)) particleEnergy[j] += 0.5*phi;
516             }
517             if (comp_energy)
518             {
519                if ((1 == HalfOrFull) && (j < numberContrib))
520                {
521                   /* Half mode -- add v to total energy */
522                   *energy += phi;
523                }
524                else
525                {
526                   /* Full mode -- add half v to total energy */
527                   *energy += 0.5*phi;
528                }
529             }
530 
531             /* contribution to process_dEdr */
532             if (comp_process_dEdr)
533             {
534                ier = KIM_API_process_dEdr(km, &dEidr, &R, &pRij, &i, &j);
535             }
536 
537             /* contribution to process_d2Edr2 */
538             if (comp_process_d2Edr2)
539             {
540                R_pairs[0] = R_pairs[1] = R;
541                Rij_pairs[0][0] = Rij_pairs[1][0] = Rij[0];
542                Rij_pairs[0][1] = Rij_pairs[1][1] = Rij[1];
543                Rij_pairs[0][2] = Rij_pairs[1][2] = Rij[2];
544                i_pairs[0] = i_pairs[1] = i;
545                j_pairs[0] = j_pairs[1] = j;
546 
547                ier = KIM_API_process_d2Edr2(km, &d2Eidr, &pR_pairs, &pRij_pairs, &pi_pairs, &pj_pairs);
548             }
549 
550             /* contribution to forces */
551             if (comp_force)
552             {
553                for (k = 0; k < DIM; ++k)
554                {
555                   force[i*DIM + k] += dEidr*Rij[k]/R; /* accumulate force on atom i */
556                   force[j*DIM + k] -= dEidr*Rij[k]/R; /* accumulate force on atom j */
557                }
558             }
559          }
560       } /* loop on jj */
561    }    /* infinite while loop (terminated by break statements above */
562 
563    /* Free temporary storage */
564    if (0 == NBC)
565    {
566       free(neighListOfCurrentAtom);
567    }
568 
569    /* everything is great */
570    ier = KIM_STATUS_OK;
571 
572    return ier;
573 }
574 
575 /* Initialization function */
_init_(void * km,char * paramfile_names,int * nmstrlen,int * numparamfiles)576 int model_driver_p_<FILL (lowercase) model driver name>_init_(void *km, char* paramfile_names, int* nmstrlen, int* numparamfiles)
577 {
578    /* KIM variables */
579    intptr_t* pkim = *((intptr_t**) km);
580    char* paramfile1name;
581    char* paramfile2name;
582    <FILL as many file name pointers as needed>
583 
584    /* Local variables */
585    FILE* fid;
586    double cutoff;
587    double <FILL parameter 1>;
588    double <FILL parameter 2>;
589    /* FILL as many parameters as needed */
590    double* model_cutoff;
591    int ier;
592    double dummy;
593    struct model_buffer* buffer;
594    char* NBCstr;
595 
596    /* set paramfile1name */
597    if (*numparamfiles != <FILL number of parameter files>)
598    {
599       ier = KIM_STATUS_FAIL;
600       KIM_API_report_error(__LINE__, __FILE__, "Incorrect number of parameter files.", ier);
601       return ier;
602    }
603    paramfile1name = paramfile_names;
604    paramfile2name = &(paramfile_names[1*(*nmstrlen)]);
605    <FILL as many file name pointers as needed>
606 
607    /* store pointer to functions in KIM object */
608    KIM_API_setm_data(pkim, &ier, 3*4,
609                      "compute", 1, &compute, 1,
610                      "reinit",  1, &reinit,  1,
611                      "destroy", 1, &destroy, 1);
612    if (KIM_STATUS_OK > ier)
613    {
614       KIM_API_report_error(__LINE__, __FILE__, "KIM_API_setm_data", ier);
615       return ier;
616    }
617 
618    /* Read in model parameters from parameter file */
619    fid = fopen(paramfile1name, "r");
620    if (fid == NULL)
621    {
622       ier = KIM_STATUS_FAIL;
623       KIM_API_report_error(__LINE__, __FILE__, "Unable to open parameter file for <FILL model driver name> parameters", ier);
624       return ier;
625    }
626 
627    ier = fscanf(fid, "%lf \n%lf \n<FILL as many parameters as needed>",
628                 &cutoff,             /* cutoff distance in angstroms */
629                 &<FILL parameter 1>,
630                 &<FILL parameter 2>,
631                 /* FILL as many parameters as needed */
632                );
633    fclose(fid);
634 
635    /* check that we read the right number of parameters */
636    if (<FILL number of parameters (including cutoff)> != ier)
637    {
638       ier = KIM_STATUS_FAIL;
639       KIM_API_report_error(__LINE__, __FILE__, "Unable to read all <FILL model driver name> parameters", ier);
640       return ier;
641    }
642 
643    /* FILL process the remaining parameter files */
644 
645    /* convert to appropriate units */
646    cutoff *= KIM_API_convert_to_act_unit(pkim, "A", "eV", "e", "K", "ps",
647                                                1.0, 0.0,  0.0, 0.0, 0.0, &ier);
648    if (KIM_STATUS_OK > ier)
649    {
650       KIM_API_report_error(__LINE__, __FILE__, "KIM_API_convert_to_act_unit", ier);
651       return ier;
652    }
653 
654    <FILL parameter 1> *= KIM_API_convert_to_act_unit(pkim, "A", "eV", "e", "K", "ps",
655                                                      <FILL exponents (5) for parameter 1>, &ier);
656    if (KIM_STATUS_OK > ier)
657    {
658       KIM_API_report_error(__LINE__, __FILE__, "KIM_API_convert_to_act_unit", ier);
659       return ier;
660    }
661 
662    <FILL parameter 2> *= KIM_API_convert_to_act_unit(pkim, "A", "eV", "e", "K", "ps",
663                                                      <FILL exponents (5) for parameter 2>, &ier);
664    if (KIM_STATUS_OK > ier)
665    {
666       KIM_API_report_error(__LINE__, __FILE__, "KIM_API_convert_to_act_unit", ier);
667       return ier;
668    }
669 
670    /* FILL as many parameters as necessary */
671 
672    /* store model cutoff in KIM object */
673    model_cutoff = (double*) KIM_API_get_data(pkim, "cutoff", &ier);
674    if (KIM_STATUS_OK > ier)
675    {
676       KIM_API_report_error(__LINE__, __FILE__, "KIM_API_get_data", ier);
677       return ier;
678    }
679    *model_cutoff = cutoff;
680 
681    /* allocate buffer */
682    buffer = (struct model_buffer*) malloc(sizeof(struct model_buffer));
683    if (NULL == buffer)
684    {
685       ier = KIM_STATUS_FAIL;
686       KIM_API_report_error(__LINE__, __FILE__, "malloc", ier);
687       return ier;
688    }
689 
690    /* setup buffer */
691    /* set value of parameters */
692    buffer->Pcutoff = *model_cutoff;
693    buffer->cutsq = (*model_cutoff)*(*model_cutoff);
694    buffer-><FILL parameter 1> = <FILL parameter 1>;
695    buffer-><FILL parameter 2> = <FILL parameter 2>;
696    /* FILL as many parameters as needed */
697 
698    /* Determine neighbor list boundary condition (NBC) */
699    NBCstr = KIM_API_get_NBC_method(pkim, &ier);
700    if (KIM_STATUS_OK > ier)
701    {
702       KIM_API_report_error(__LINE__, __FILE__, "KIM_API_get_NBC_method", ier);
703       return ier;
704    }
705    if (!strcmp("CLUSTER",NBCstr))
706    {
707       buffer->NBC = 0;
708    }
709    else if ((!strcmp("MI_OPBC_H",NBCstr)) || (!strcmp("MI_OPBC_F",NBCstr)))
710    {
711       buffer->NBC = 1;
712    }
713    else if ((!strcmp("NEIGH_PURE_H",NBCstr)) || (!strcmp("NEIGH_PURE_F",NBCstr)))
714    {
715       buffer->NBC = 2;
716    }
717    else if (!strcmp("NEIGH_RVEC_F",NBCstr))
718    {
719       buffer->NBC = 3;
720    }
721    else
722    {
723       ier = KIM_STATUS_FAIL;
724       KIM_API_report_error(__LINE__, __FILE__, "Unknown NBC method", ier);
725       return ier;
726    }
727    free(NBCstr); /* don't forget to release the memory... */
728 
729    /* Determine if Half or Full neighbor lists are being used */
730    /*****************************
731     * HalfOrFull = 1 -- Half
732     *            = 2 -- Full
733     *****************************/
734    if (KIM_API_is_half_neighbors(pkim, &ier))
735    {
736       buffer->HalfOrFull = 1;
737    }
738    else
739    {
740       buffer->HalfOrFull = 2;
741    }
742 
743    /* determine neighbor list handling mode */
744    if (buffer->NBC != 0)
745    {
746       /*****************************
747        * IterOrLoca = 1 -- Iterator
748        *            = 2 -- Locator
749        *****************************/
750       buffer->IterOrLoca = KIM_API_get_neigh_mode(pkim, &ier);
751       if (KIM_STATUS_OK > ier)
752       {
753          KIM_API_report_error(__LINE__, __FILE__, "KIM_API_get_neigh_mode", ier);
754          return ier;
755       }
756       if ((buffer->IterOrLoca != 1) && (buffer->IterOrLoca != 2))
757       {
758          ier = KIM_STATUS_FAIL;
759          KIM_API_report_error(__LINE__, __FILE__, "Unsupported IterOrLoca mode", ier);
760          return ier;
761       }
762    }
763    else
764    {
765       buffer->IterOrLoca = 2;   /* for CLUSTER NBC */
766    }
767 
768    buffer->model_index_shift = KIM_API_get_model_index_shift(pkim);
769 
770    KIM_API_getm_index(pkim, &ier, 12*3,
771                       "cutoff",                      &(buffer->cutoff_ind),                      1,
772                       "numberOfParticles",           &(buffer->numberOfParticles_ind),           1,
773                       "particleTypes",               &(buffer->particleTypes_ind),               1,
774                       "numberContributingParticles", &(buffer->numberContributingParticles_ind), 1,
775                       "coordinates",                 &(buffer->coordinates_ind),                 1,
776                       "get_neigh",                   &(buffer->get_neigh_ind),                   1,
777                       "boxSideLengths",              &(buffer->boxSideLengths_ind),              1,
778                       "energy",                      &(buffer->energy_ind),                      1,
779                       "forces",                      &(buffer->forces_ind),                      1,
780                       "particleEnergy",              &(buffer->particleEnergy_ind),              1,
781                       "process_dEdr",                &(buffer->process_dEdr_ind),                1,
782                       "process_d2Edr2",              &(buffer->process_d2Edr2_ind),              1
783                      );
784    if (KIM_STATUS_OK > ier)
785    {
786       KIM_API_report_error(__LINE__, __FILE__, "KIM_API_getm_index", ier);
787       return ier;
788    }
789    /* end setup buffer */
790 
791    /* store in model buffer */
792    KIM_API_set_model_buffer(pkim, (void*) buffer, &ier);
793    if (KIM_STATUS_OK > ier)
794    {
795       KIM_API_report_error(__LINE__, __FILE__, "KIM_API_set_model_buffer", ier);
796       return ier;
797    }
798 
799    /* set pointers to parameters in KIM object */
800    KIM_API_setm_data(pkim, &ier, 6*4,
801                      "PARAM_FREE_cutoff",  1, &(buffer->Pcutoff), 1,
802                      "PARAM_FIXED_cutsq",  1, &(buffer->cutsq),   1,
803                      "PARAM_<FILL FREE or FIXED>_<FILL parameter 1>", 1, &(buffer-><FILL parameter 1>), 1,
804                      "PARAM_<FILL FREE or FIXED>_<FILL parameter 2>", 1, &(buffer-><FILL parameter 2>), 1,
805                      /* FILL as many parameters as needed */
806                     );
807    if (KIM_STATUS_OK > ier)
808    {
809       KIM_API_report_error(__LINE__, __FILE__, "KIM_API_setm_data", ier);
810       return ier;
811    }
812 
813    return KIM_STATUS_OK;
814 }
815 
816 /* Reinitialization function */
reinit(void * km)817 static int reinit(void *km)
818 {
819    /* Local variables */
820    intptr_t* pkim = *((intptr_t**) km);
821    int ier;
822    double *cutoff;
823    double dummy;
824    struct model_buffer* buffer;
825 
826    /* get buffer from KIM object */
827    buffer = (struct model_buffer*) KIM_API_get_model_buffer(pkim, &ier);
828    if (KIM_STATUS_OK > ier)
829    {
830       KIM_API_report_error(__LINE__, __FILE__, "KIM_API_get_model_buffer", ier);
831       return ier;
832    }
833 
834    /* set new values in KIM object     */
835    /*                                  */
836    /* store model cutoff in KIM object */
837    cutoff = KIM_API_get_data(pkim, "cutoff", &ier);
838    if (KIM_STATUS_OK > ier)
839    {
840       KIM_API_report_error(__LINE__, __FILE__, "KIM_API_get_data", ier);
841       return ier;
842    }
843    *cutoff = buffer->Pcutoff;
844 
845    /* set value of parameter cutsq */
846    buffer->cutsq = (*cutoff)*(*cutoff);
847 
848    /* FILL: recompute any other FIXED parameters whose values depend on FREE parameters */
849 
850    ier = KIM_STATUS_OK;
851    return ier;
852 }
853 
854 /* destroy function */
destroy(void * km)855 static int destroy(void *km)
856 {
857    /* Local variables */
858    intptr_t* pkim = *((intptr_t**) km);
859    struct model_buffer* buffer;
860    int ier;
861 
862    /* get model buffer from KIM object */
863    buffer = (struct model_buffer*) KIM_API_get_model_buffer(pkim, &ier);
864    if (KIM_STATUS_OK > ier)
865    {
866       KIM_API_report_error(__LINE__, __FILE__, "KIM_API_get_model_buffer", ier);
867       return ier;
868    }
869 
870    /* destroy the buffer */
871    free(buffer);
872 
873    ier = KIM_STATUS_OK;
874    return ier;
875 }
876