KIM API V2
ex_test_Ar_fcc_cluster.c
Go to the documentation of this file.
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) 2013--2018, Regents of the University of Minnesota.
24  * All rights reserved.
25  *
26  * Contributors:
27  * Ryan S. Elliott
28  * Stephen M. Whalen
29  *
30  */
31 
32 
33 #include <stdlib.h>
34 #include <stdio.h>
35 #include <math.h>
36 #include "KIM_SimulatorHeaders.h"
37 
38 #define TRUE 1
39 #define FALSE 0
40 
41 #define NAMESTRLEN 128
42 
43 #define FCCSPACING 5.260
44 #define DIM 3
45 #define NCELLSPERSIDE 2
46 #define NCLUSTERPARTS (4*(NCELLSPERSIDE*NCELLSPERSIDE*NCELLSPERSIDE) + \
47  6*(NCELLSPERSIDE*NCELLSPERSIDE) \
48  + 3*(NCELLSPERSIDE) + 1)
49 
50 
51 #define MY_ERROR(message) \
52  { \
53  printf("* Error : \"%s\" %d:%s\n", message, \
54  __LINE__, __FILE__); \
55  exit(1); \
56  }
57 
58 #define MY_WARNING(message) \
59  { \
60  printf("* Error : \"%s\" %d:%s\n", message, \
61  __LINE__, __FILE__); \
62  }
63 
64 /* Define neighborlist structure */
65 typedef struct
66 {
67  double cutoff;
69  int* NNeighbors;
70  int* neighborList;
71 } NeighList;
72 
73 /* Define prototypes */
74 void fcc_cluster_neighborlist(int half, int numberOfParticles, double* coords,
75  double cutoff, NeighList* nl);
76 
78  void const * const dataObject,
79  int const numberOfNeighborLists, double const * const cutoffs,
80  int const neighborListIndex,
81  int const particleNumber, int * const numberOfNeighbors,
82  int const ** const neighborsOfParticle);
83 
84 void create_FCC_cluster(double FCCspacing, int nCellsPerSide, double *coords);
85 
86 
87 /* Main program */
88 int main()
89 {
90  /* Local variable declarations */
91  double const MinSpacing = 0.8*FCCSPACING;
92  double const MaxSpacing = 1.2*FCCSPACING;
93  double const SpacingIncr = 0.025*FCCSPACING;
94  double CurrentSpacing;
95  double cutpad = 0.75; /* Angstroms */
96  double force_norm;
97  int i;
98  int error;
99 
100 
101  /* KIM variable declarations */
102  KIM_Model * model;
103  /* model inputs */
104  int numberOfParticles_cluster = NCLUSTERPARTS;
105  int speciesIsSupported;
106  int particleSpecies_cluster_model[NCLUSTERPARTS];
107  int particleContributing_cluster_model[NCLUSTERPARTS];
108  double coords_cluster[NCLUSTERPARTS][DIM];
109  NeighList nl_cluster_model;
110  /* model outputs */
111  int number_of_neighbor_lists_cluster_model;
112  double influence_distance_cluster_model;
113  double const * cutoff_cluster_model;
114  double energy_cluster_model;
115  double forces_cluster[NCLUSTERPARTS*DIM];
116 
117  char modelname[NAMESTRLEN];
118  int requestedUnitsAccepted;
119  int modelArCode;
120  int numberOfComputeArgumentNames;
121  KIM_ComputeArguments * computeArguments;
122  KIM_ComputeArgumentName computeArgumentName;
123  KIM_SupportStatus supportStatus;
124  int numberOfComputeCallbackNames;
125  KIM_ComputeCallbackName computeCallbackName;
126 
127  /* Get KIM Model names */
128  printf("Please enter valid KIM Model name: \n");
129  error = scanf("%s", modelname);
130  if (1 != error)
131  {
132  MY_ERROR("Unable to read model name");
133  }
134 
135  /* initialize the model */
142  modelname,
143  &requestedUnitsAccepted,
144  &model);
145  if (error) MY_ERROR("KIM_create_model_interface()");
146 
147  /* Check for compatibility with the model */
148  if (!requestedUnitsAccepted) MY_ERROR("Must adapt to model units");
149 
150  /* check species */
152  model, KIM_SPECIES_NAME_Ar, &speciesIsSupported, &modelArCode);
153  if ((error) || (!speciesIsSupported))
154  {
155  MY_ERROR("Species Ar not supported");
156  }
157 
158  error = KIM_Model_ComputeArgumentsCreate(model, &computeArguments);
159  if (error)
160  {
161  MY_ERROR("KIM_Model_ComputeArgumentsCreate");
162  }
163 
164  /* check arguments */
166  &numberOfComputeArgumentNames);
167  for (i=0; i<numberOfComputeArgumentNames; ++i)
168  {
170  i, &computeArgumentName);
171  if (error) MY_ERROR("can't get argument name");
173  computeArgumentName,
174  &supportStatus);
175  if (error) MY_ERROR("can't get argument supportStatus");
176 
177  /* can only handle energy and force as a required arg */
179  {
181  computeArgumentName,
184  computeArgumentName,
186  {
187  MY_ERROR("unsupported required argument");
188  }
189  }
190 
191  /* must have energy and forces */
193  computeArgumentName,
195  ||
197  computeArgumentName,
199  {
200  if (! ((KIM_SupportStatus_Equal(supportStatus,
202  ||
203  (KIM_SupportStatus_Equal(supportStatus,
205  {
206  MY_ERROR("energy or forces not available");
207  }
208  }
209  }
210 
211  /* check call backs */
213  &numberOfComputeCallbackNames);
214  for (i=0; i<numberOfComputeCallbackNames; ++i)
215  {
217  i, &computeCallbackName);
218  if (error) MY_ERROR("can't get call back name");
220  computeCallbackName,
221  &supportStatus);
222  if (error) MY_ERROR("can't get call back supportStatus");
223 
224  /* cannot handle any "required" call backs */
226  {
227  MY_ERROR("unsupported required call back");
228  }
229  }
230 
231  /* We're compatible with the model. Let's do it. */
232 
233  error =
235  computeArguments,
237  &numberOfParticles_cluster)
238  ||
240  computeArguments,
242  particleSpecies_cluster_model)
243  ||
245  computeArguments,
247  particleContributing_cluster_model)
248  ||
250  computeArguments, KIM_COMPUTE_ARGUMENT_NAME_coordinates,
251  (double*) &coords_cluster)
252  ||
254  computeArguments, KIM_COMPUTE_ARGUMENT_NAME_partialEnergy,
255  &energy_cluster_model)
256  ||
258  computeArguments, KIM_COMPUTE_ARGUMENT_NAME_partialForces,
259  (double*) &forces_cluster);
260  if (error) MY_ERROR("KIM_setm_data");
262  computeArguments,
265  (func *) &get_cluster_neigh,
266  &nl_cluster_model);
267 
268  KIM_Model_GetInfluenceDistance(model, &influence_distance_cluster_model);
270  &number_of_neighbor_lists_cluster_model,
271  &cutoff_cluster_model,
272  NULL, /* ignoring hints */
273  NULL); /* ignoring hints */
274  if (number_of_neighbor_lists_cluster_model != 1)
275  MY_ERROR("too many neighbor lists");
276 
277  /* setup particleSpecies */
279  model,
281  &speciesIsSupported,
282  &(particleSpecies_cluster_model[0]));
283  if (error) MY_ERROR("KIM_get_species_code");
284  for (i = 1; i < NCLUSTERPARTS; ++i)
285  particleSpecies_cluster_model[i] = particleSpecies_cluster_model[0];
286 
287  /* setup particleContributing */
288  for (i = 0; i < NCLUSTERPARTS; ++i)
289  {
290  particleContributing_cluster_model[i] = 1; /* all particles contribute */
291  }
292 
293  /* setup neighbor lists */
294  /* allocate memory for list */
295  nl_cluster_model.numberOfParticles = NCLUSTERPARTS;
296  nl_cluster_model.NNeighbors = (int*) malloc(NCLUSTERPARTS*sizeof(int));
297  if (NULL==nl_cluster_model.NNeighbors) MY_ERROR("malloc unsuccessful");
298 
299  nl_cluster_model.neighborList = (int*) malloc(NCLUSTERPARTS*NCLUSTERPARTS*sizeof(int));
300  if (NULL==nl_cluster_model.neighborList) MY_ERROR("malloc unsuccessful");
301 
302  /* ready to compute */
303  printf("This is Test : ex_test_Ar_fcc_cluster\n");
304  printf("--------------------------------------------------------------------------------\n");
305  printf("Results for KIM Model : %s\n", modelname);
306 
307  printf("%20s, %20s, %20s\n", "Energy", "Force Norm", "Lattice Spacing");
308  for (CurrentSpacing = MinSpacing; CurrentSpacing < MaxSpacing; CurrentSpacing += SpacingIncr)
309  {
310  /* update coordinates for cluster */
311  create_FCC_cluster(CurrentSpacing, NCELLSPERSIDE, &(coords_cluster[0][0]));
312  /* compute neighbor lists */
313  fcc_cluster_neighborlist(0, NCLUSTERPARTS, &(coords_cluster[0][0]),
314  (*cutoff_cluster_model + cutpad), &nl_cluster_model);
315 
316  /* call compute functions */
317  error = KIM_Model_Compute(model, computeArguments);
318  if (error) MY_ERROR("KIM_model_compute");
319 
320  /* compute force norm */
321  force_norm = 0.0;
322  for (i=0; i < DIM*numberOfParticles_cluster; ++i)
323  {
324  force_norm += forces_cluster[i]*forces_cluster[i];
325  }
326  force_norm = sqrt(force_norm);
327 
328  /* print the results */
329  printf("%20.10e, %20.10e, %20.10e\n",
330  energy_cluster_model,
331  force_norm,
332  CurrentSpacing);
333  }
334 
335  error = KIM_Model_ComputeArgumentsDestroy(model, &computeArguments);
336  if (error)
337  {
338  MY_ERROR("Unable to destroy compute arguments");
339  }
340 
341  /* free memory of neighbor lists */
342  free(nl_cluster_model.NNeighbors);
343  free(nl_cluster_model.neighborList);
344 
345  /* free pkim objects */
346  KIM_Model_Destroy(&model);
347 
348  /* everything is great */
349  return 0;
350 }
351 
352 void create_FCC_cluster(double FCCspacing, int nCellsPerSide, double *coords)
353 {
354  /* local variables */
355  double FCCshifts[4][DIM];
356  double latVec[DIM];
357  int a;
358  int i;
359  int j;
360  int k;
361  int m;
362  int n;
363 
364  /* create a cubic FCC cluster of parts */
365  FCCshifts[0][0] = 0.0; FCCshifts[0][1] = 0.0; FCCshifts[0][2] = 0.0;
366  FCCshifts[1][0] = 0.5*FCCspacing; FCCshifts[1][1] = 0.5*FCCspacing; FCCshifts[1][2] = 0.0;
367  FCCshifts[2][0] = 0.5*FCCspacing; FCCshifts[2][1] = 0.0; FCCshifts[2][2] = 0.5*FCCspacing;
368  FCCshifts[3][0] = 0.0; FCCshifts[3][1] = 0.5*FCCspacing; FCCshifts[3][2] = 0.5*FCCspacing;
369 
370  a = 0;
371  for (i = 0; i < nCellsPerSide; ++i)
372  {
373  latVec[0] = ((double) i)*FCCspacing;
374  for (j = 0; j < nCellsPerSide; ++j)
375  {
376  latVec[1] = ((double) j)*FCCspacing;
377  for (k = 0; k < nCellsPerSide; ++k)
378  {
379  latVec[2] = ((double) k)*FCCspacing;
380  for (m = 0; m < 4; ++m)
381  {
382  for (n = 0; n < DIM; ++n)
383  {
384  coords[a*DIM + n] = latVec[n] + FCCshifts[m][n];
385  }
386  a++;
387  }
388  }
389  /* add in the remaining three faces */
390  /* pos-x face */
391  latVec[0] = NCELLSPERSIDE*FCCspacing;
392  latVec[1] = ((double) i)*FCCspacing;
393  latVec[2] = ((double) j)*FCCspacing;
394  for (n = 0; n < DIM; ++n)
395  {
396  coords[a*DIM + n] = latVec[n];
397  }
398  a++;
399  for (n = 0; n < DIM; ++n)
400  {
401  coords[a*DIM + n] = latVec[n] + FCCshifts[3][n];
402  }
403  a++;
404  /* pos-y face */
405  latVec[0] = ((double) i)*FCCspacing;
406  latVec[1] = NCELLSPERSIDE*FCCspacing;
407  latVec[2] = ((double) j)*FCCspacing;
408  for (n = 0; n < DIM; ++n)
409  {
410  coords[a*DIM + n] = latVec[n];
411  }
412  a++;
413  for (n = 0; n < DIM; ++n)
414  {
415  coords[a*DIM + n] = latVec[n] + FCCshifts[2][n];
416  }
417  a++;
418  /* pos-z face */
419  latVec[0] = ((double) i)*FCCspacing;
420  latVec[1] = ((double) j)*FCCspacing;
421  latVec[2] = NCELLSPERSIDE*FCCspacing;
422  for (n = 0; n < DIM; ++n)
423  {
424  coords[a*DIM + n] = latVec[n];
425  }
426  a++;
427  for (n = 0; n < DIM; ++n)
428  {
429  coords[a*DIM + n] = latVec[n] + FCCshifts[1][n];
430  }
431  a++;
432  }
433  /* add in the remaining three edges */
434  latVec[0] = ((double) i)*FCCspacing;
435  latVec[1] = NCELLSPERSIDE*FCCspacing;
436  latVec[2] = NCELLSPERSIDE*FCCspacing;
437  for (n = 0; n < DIM; ++n)
438  {
439  coords[a*DIM + n] = latVec[n];
440  }
441  a++;
442  latVec[0] = NCELLSPERSIDE*FCCspacing;
443  latVec[1] = ((double) i)*FCCspacing;
444  latVec[2] = NCELLSPERSIDE*FCCspacing;
445  for (n = 0; n < DIM; ++n)
446  {
447  coords[a*DIM + n] = latVec[n];
448  }
449  a++;
450  latVec[0] = NCELLSPERSIDE*FCCspacing;
451  latVec[1] = NCELLSPERSIDE*FCCspacing;
452  latVec[2] = ((double) i)*FCCspacing;
453  for (n = 0; n < DIM; ++n)
454  {
455  coords[a*DIM + n] = latVec[n];
456  }
457  a++;
458  }
459  /* add in the remaining corner */
460  for (n = 0; n < DIM; ++n)
461  {
462  coords[a*DIM + n] = NCELLSPERSIDE*FCCspacing;
463  }
464  a++;
465 
466  return;
467 }
468 
469 
470 void fcc_cluster_neighborlist(int half, int numberOfParticles, double* coords,
471  double cutoff, NeighList* nl)
472 {
473  /* local variables */
474  int i;
475  int j;
476  int k;
477  int a;
478 
479  double dx[DIM];
480  double r2;
481  double cutoff2;
482 
483  nl->cutoff = cutoff;
484 
485  cutoff2 = cutoff*cutoff;
486 
487  for (i = 0; i < numberOfParticles; ++i)
488  {
489  a = 0;
490  for (j = 0; j < numberOfParticles; ++j)
491  {
492  r2 = 0.0;
493  for (k = 0; k < DIM; ++k)
494  {
495  dx[k] = coords[j*DIM + k] - coords[i*DIM + k];
496  r2 += dx[k]*dx[k];
497  }
498 
499  if (r2 < cutoff2)
500  {
501  if ((half && i < j) || (!half && i != j))
502  {
503  /* part j is a neighbor of part i */
504  (*nl).neighborList[i*NCLUSTERPARTS + a] = j;
505  a++;
506  }
507  }
508  }
509  /* part i has `a' neighbors */
510  (*nl).NNeighbors[i] = a;
511  }
512 
513  return;
514 }
515 
516 int get_cluster_neigh(void const * const dataObject,
517  int const numberOfNeighborLists,
518  double const * const cutoffs,
519  int const neighborListIndex,
520  int const particleNumber, int * const numberOfNeighbors,
521  int const ** const neighborsOfParticle)
522 {
523  /* local variables */
524  int error = TRUE;
525  NeighList* nl = (NeighList*) dataObject;
527 
528  if ((numberOfNeighborLists != 1) || (cutoffs[0] > nl->cutoff)) return error;
529 
530  if (neighborListIndex != 0) return error;
531 
532  /* initialize numNeigh */
533  *numberOfNeighbors = 0;
534 
535  if ((particleNumber >= numberOfParticles) || (particleNumber < 0)) /* invalid id */
536  {
537  MY_WARNING("Invalid part ID in get_cluster_neigh");
538  return TRUE;
539  }
540 
541  /* set the returned number of neighbors for the returned part */
542  *numberOfNeighbors = (*nl).NNeighbors[particleNumber];
543 
544  /* set the location for the returned neighbor list */
545  *neighborsOfParticle = &((*nl).neighborList[(particleNumber)*numberOfParticles]);
546 
547  return FALSE;
548 }
#define DIM
struct KIM_ComputeArguments KIM_ComputeArguments
int get_cluster_neigh(void *kimmdl, int *mode, int *request, int *part, int *numnei, int **nei1part, double **Rij)
void KIM_Model_GetInfluenceDistance(KIM_Model const *const model, double *const influenceDistance)
KIM_SupportStatus const KIM_SUPPORT_STATUS_optional
int KIM_ComputeArguments_GetCallbackSupportStatus(KIM_ComputeArguments const *const computeArguments, KIM_ComputeCallbackName const computeCallbackName, KIM_SupportStatus *const supportStatus)
void KIM_COMPUTE_ARGUMENT_NAME_GetNumberOfComputeArgumentNames(int *const numberOfComputeArgumentNames)
KIM_ComputeArgumentName const KIM_COMPUTE_ARGUMENT_NAME_coordinates
void KIM_COMPUTE_CALLBACK_NAME_GetNumberOfComputeCallbackNames(int *const numberOfComputeCallbackNames)
void KIM_Model_Destroy(KIM_Model **const model)
#define MY_ERROR(message)
KIM_ChargeUnit const KIM_CHARGE_UNIT_e
KIM_LanguageName const KIM_LANGUAGE_NAME_c
int KIM_ComputeArguments_SetCallbackPointer(KIM_ComputeArguments *const computeArguments, KIM_ComputeCallbackName const computeCallbackName, KIM_LanguageName const languageName, func *const fptr, void const *const dataObject)
KIM_ComputeArgumentName const KIM_COMPUTE_ARGUMENT_NAME_partialForces
KIM_ComputeArgumentName const KIM_COMPUTE_ARGUMENT_NAME_particleSpeciesCodes
#define FCCSPACING
int KIM_Model_GetSpeciesSupportAndCode(KIM_Model const *const model, KIM_SpeciesName const speciesName, int *const speciesIsSupported, int *const code)
#define NCELLSPERSIDE
int KIM_ComputeArguments_SetArgumentPointerInteger(KIM_ComputeArguments *const computeArguments, KIM_ComputeArgumentName const computeArgumentName, int const *const ptr)
KIM_SpeciesName const KIM_SPECIES_NAME_Ar
int KIM_Model_ComputeArgumentsCreate(KIM_Model const *const model, KIM_ComputeArguments **const computeArguments)
#define FALSE
void() func()
Definition: KIM_func.h:39
#define NAMESTRLEN
LengthUnit const m
int KIM_Model_Create(KIM_Numbering const numbering, KIM_LengthUnit const requestedLengthUnit, KIM_EnergyUnit const requestedEnergyUnit, KIM_ChargeUnit const requestedChargeUnit, KIM_TemperatureUnit const requestedTemperatureUnit, KIM_TimeUnit const requestedTimeUnit, char const *const modelName, int *const requestedUnitsAccepted, KIM_Model **const model)
KIM_ComputeArgumentName const KIM_COMPUTE_ARGUMENT_NAME_numberOfParticles
int KIM_ComputeArgumentName_NotEqual(KIM_ComputeArgumentName const left, KIM_ComputeArgumentName const right)
KIM_ComputeCallbackName const KIM_COMPUTE_CALLBACK_NAME_GetNeighborList
#define TRUE
KIM_TimeUnit const KIM_TIME_UNIT_ps
int KIM_ComputeArguments_SetArgumentPointerDouble(KIM_ComputeArguments *const computeArguments, KIM_ComputeArgumentName const computeArgumentName, double const *const ptr)
int KIM_COMPUTE_CALLBACK_NAME_GetComputeCallbackName(int const index, KIM_ComputeCallbackName *const computeCallbackName)
void create_FCC_cluster(double FCCspacing, int nCellsPerSide, double *coords)
int KIM_Model_ComputeArgumentsDestroy(KIM_Model const *const model, KIM_ComputeArguments **const computeArguments)
int KIM_SupportStatus_Equal(KIM_SupportStatus const left, KIM_SupportStatus const right)
KIM_EnergyUnit const KIM_ENERGY_UNIT_eV
#define MY_WARNING(message)
#define NCLUSTERPARTS
ComputeArgumentName const numberOfParticles
int KIM_Model_Compute(KIM_Model const *const model, KIM_ComputeArguments const *const computeArguments)
int KIM_COMPUTE_ARGUMENT_NAME_GetComputeArgumentName(int const index, KIM_ComputeArgumentName *const computeArgumentName)
KIM_ComputeArgumentName const KIM_COMPUTE_ARGUMENT_NAME_partialEnergy
void KIM_Model_GetNeighborListPointers(KIM_Model const *const model, int *const numberOfNeighborLists, double const **const cutoffs, int const **const paddingNeighborHints, int const **const halfListHints)
KIM_SupportStatus const KIM_SUPPORT_STATUS_required
int main()
int KIM_ComputeArguments_GetArgumentSupportStatus(KIM_ComputeArguments const *const computeArguments, KIM_ComputeArgumentName const computeArgumentName, KIM_SupportStatus *const supportStatus)
struct KIM_Model KIM_Model
Definition: KIM_Model.h:104
void fcc_cluster_neighborlist(int half, int numberOfParticles, double *coords, double cutoff, NeighList *nl)
KIM_LengthUnit const KIM_LENGTH_UNIT_A
LogVerbosity const error
int KIM_ComputeArgumentName_Equal(KIM_ComputeArgumentName const left, KIM_ComputeArgumentName const right)
KIM_TemperatureUnit const KIM_TEMPERATURE_UNIT_K
KIM_Numbering const KIM_NUMBERING_zeroBased
KIM_ComputeArgumentName const KIM_COMPUTE_ARGUMENT_NAME_particleContributing