Added code to enable generation of multilayers and statistics collect…

…ion on a layerwise basis

hdr/create.hpp

@@ -98,6 +98,11 @@ namespace cs{
 	extern double interfacial_roughness_mean_seed_height; // Angstroms
 	extern double interfacial_roughness_seed_height_max; // Angstroms
 
+   // Variables for multilayer system
+   extern bool multilayers;
+   extern int num_multilayers;
+   extern bool multilayer_height_category; // enable height categorization by multilayer number
+
 	class unit_cell_atom_t {
 	public:
 		double x; /// atom x-coordinate
@@ -438,6 +443,7 @@ int sort_atoms_by_grain(std::vector<cs::catom_t> &);
 int clear_atoms(std::vector<cs::catom_t> &);
 
 void roughness(std::vector<cs::catom_t> &);
+void generate_multilayers(std::vector<cs::catom_t> & catom_array);
 
   // unit cell initialisation function
   void unit_cell_set(cs::unit_cell_t &);

makefile

@@ -65,6 +65,7 @@ obj/create/cs_create_neighbour_list2.o \
 obj/create/cs_particle_shapes.o \
 obj/create/cs_set_atom_vars2.o \
 obj/create/cs_voronoi2.o \
+obj/create/multilayers.o \
 obj/data/atoms.o \
 obj/data/category.o \
 obj/data/cells.o \

src/create/create_system2.cpp

@@ -98,6 +98,11 @@ namespace cs{
 	std::string unit_cell_file="";
 	bool fill_core_shell=true;
 
+   // Variables for multilayer system
+   bool multilayers = false;
+   bool multilayer_height_category = false; // enable height categorization by multilayer number
+   int num_multilayers = 1;
+
 	// Variables for interfacial roughness control
 	bool interfacial_roughness=false;
 	bool interfacial_roughness_local_height_field=false;

src/create/cs_create_crystal_structure2.cpp

@@ -172,6 +172,9 @@ int create_crystal_structure(std::vector<cs::catom_t> & catom_array){
 		// Check for interfacial roughness and call custom material assignment routine
 		if(cs::interfacial_roughness==true) cs::roughness(catom_array);
 
+      // Check for multilayer system and if required generate multilayers
+      else if(cs::multilayers) cs::generate_multilayers(catom_array);
+
 		// Otherwise perform normal assignement of materials
 		else{
 

src/create/multilayers.cpp

@@ -0,0 +1,79 @@
+//-----------------------------------------------------------------------------
+//
+// This source file is part of the VAMPIRE open source package under the
+// GNU GPL (version 2) licence (see licence file for details).
+//
+// (c) R F L Evans 2014. All rights reserved.
+//
+//-----------------------------------------------------------------------------
+
+// C++ standard library headers
+
+// Vampire headers
+#include "create.hpp"
+#include "material.hpp"
+
+namespace cs{
+
+   //--------------------------------------------------------------------------
+   // Function to create a multilayer system replicating a standard stack 
+   // structure cs::num_multilayers times
+   //
+   //   |-----------------------|  1.0          |-----------------------| L2 M2
+   //   |      Material 2       |               |-----------------------|
+   //   |-----------------------|  0.6   x2     |                       | L2 M2
+   //   |                       |        -->    |-----------------------| L1 M1
+   //   |      Material 1       |               |-----------------------|
+   //   |                       |               |                       | L1 M2 
+   //   |-----------------------|  0.0          |-----------------------|
+   //
+   // The multilayers code divides the total system height into n multiples
+   // of the defined material heights.
+   //--------------------------------------------------------------------------
+   void generate_multilayers(std::vector<cs::catom_t>& catom_array){
+
+      // Load number of multilayer repeats to temporary constant
+      const int num_layers = cs::num_multilayers;
+
+      // Determine fractional system height for each layer
+      const double fractional_system_height = cs::system_dimensions[2]/double(num_layers);
+
+      // Unroll min, max and fill for performance
+      std::vector<double> mat_min(mp::num_materials);
+      std::vector<double> mat_max(mp::num_materials);
+      std::vector<bool> mat_fill(mp::num_materials);
+
+      for(int mat=0;mat<mp::num_materials;mat++){
+         mat_min[mat]=mp::material[mat].min;
+         mat_max[mat]=mp::material[mat].max;
+         // alloys generally are not defined by height, and so have max = 0.0
+         if(mat_max[mat]<0.0000001) mat_max[mat]=-0.1;
+         mat_fill[mat]=mp::material[mat].fill;
+      }
+
+      // Assign materials to generated atoms
+      for(unsigned int atom=0;atom<catom_array.size();atom++){
+         // loop over multilayers
+         for(int multi=0; multi < num_layers; ++multi){
+            const double multilayer_num = double(multi);
+            for(int mat=0;mat<mp::num_materials;mat++){
+               // determine mimimum and maximum heigh for this layer
+               const double mat_min_z = (mat_min[mat] + multilayer_num)*fractional_system_height;
+               const double mat_max_z = (mat_max[mat] + multilayer_num)*fractional_system_height;
+               const double cz=catom_array[atom].z;
+               // if in range then allocate to material
+               if((cz>=mat_min_z) && (cz<mat_max_z) && (mat_fill[mat]==false)){
+                  catom_array[atom].material=mat;
+                  catom_array[atom].include=true;
+                  // Optionally recategorize heigh magnetization by layer in multilayer
+                  if(cs::multilayer_height_category) catom_array[atom].lh_category = multi;
+               }
+            }
+         }
+      }
+
+      return;
+
+   }
+
+} // end of namespace cs

src/utility/vio.cpp

@@ -1176,6 +1176,42 @@ int match_create(string const word, string const value, string const unit, int c
       return EXIT_SUCCESS;
    }
    //--------------------------------------------------------------------
+   test="multilayers";
+   if(word==test){
+      int nmul=atoi(value.c_str());
+      // Test for valid range
+      check_for_valid_int(nmul, word, line, prefix, 1, 100,"input","1 - 100, specifying the number of multilayers to be generated");
+      cs::multilayers = true;
+      cs::num_multilayers = nmul;
+      return EXIT_SUCCESS;
+   }
+   //--------------------------------------------------------------------
+   test="height-categorization";
+   if(word==test){
+      // Test for different options
+      test="default";
+      if(value==test){
+         // do nothing
+         return EXIT_SUCCESS;
+      }
+      test="multilayers";
+      if(value==test){
+         cs::multilayer_height_category = true;
+         return EXIT_SUCCESS;
+      }
+      else{
+         terminaltextcolor(RED);
+         std::cerr << "Error - value for \'create:" << word << "\' must be one of:" << std::endl;
+         std::cerr << "\t\"default\"" << std::endl;
+         std::cerr << "\t\"multilayers\"" << std::endl;
+         zlog << zTs() << "Error - value for \'create:" << word << "\' must be one of:" << std::endl;
+         zlog << zTs() << "\t\"default\"" << std::endl;
+         zlog << zTs() << "\t\"multilayers\"" << std::endl;
+         terminaltextcolor(WHITE);
+         err::vexit();
+      }
+   }
+   //--------------------------------------------------------------------
    // keyword not found
    //--------------------------------------------------------------------
    else{