oneoff_optimize_1d_argprocess.c++

Go to the documentation of this file.

/*
 *oneoff_1d_argprocess: process args for oneoff_1d

    Copyright (C) 2010 Joseph Pingenot

    This program is free software: you can redistribute it and/or modify
    it under the terms of the GNU Affero General Public License as published by
    the Free Software Foundation, either version 3 of the License, or
    (at your option) any later version.

    This program is distributed in the hope that it will be useful,
    but WITHOUT ANY WARRANTY; without even the implied warranty of
    MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
    GNU Affero General Public License for more details.

    You should have received a copy of the GNU Affero General Public License
    along with this program.  If not, see <http://www.gnu.org/licenses/>.

*/

/*The following will be stripped out in the MPI-enabled file*/
#undef HAVE_MPI

#include "oneoff_optimize_1d_argprocess.h++"
#include <math.h>
#include <stdio.h>
#include <glib.h>
#include <gio/gio.h>
#include <gio/gunixinputstream.h>
#include <stdlib.h>
#include "kdotp/makestruct_argprocess.h++"
#include <libmodelxx/function_specification.h++>
#include "config.h"

using namespace kdotp::libmodelxx::postprocessing;

extern "C" {
  void program_print_version(FILE* f);
  void libkdotp_print_version(FILE* f);
  void libmodelxx_print_version(FILE* f);
  void liblapackwrap_print_version(FILE* f);
  void libpostproc_print_version(FILE* f);
  void libpostproccalx_print_version(FILE* f);
}

void print_version_info(void) {
  program_print_version(stdout);
  libmodelxx_print_version(stdout);
  libkdotp_print_version(stdout);
  liblapackwrap_print_version(stdout);
  libpostproc_print_version(stdout);
  libpostproc_print_version(stdout);
  libpostproccalx_print_version(stdout);
  printf("autotools info:\n\tProject: %s\n\tVersion: %s\n", PACKAGE_NAME, VERSION);
}

struct template_parameter_range* new_template_parameter_range(char* name, double from, double to, double step) {
  struct template_parameter_range *tpr = (struct template_parameter_range*)malloc(sizeof(template_parameter_range));
  if(tpr == NULL) return NULL;
  tpr->symbol = name;
  tpr->from = from;
  tpr->to = to;
  tpr->stepsize = step;
  return tpr;
}

/**Returns a pointer to a malloc'ed double if successful,
 *  else returns NULL
 */
double* get_double_for_string(const char* s) {
  double v;
  int n = sscanf(s, "%lg", &v);
  if(n != 1) return NULL;
  double* pv = (double*)malloc(sizeof(double));
  if(pv == NULL) return NULL;
  *pv = v;
  return pv;
}
unsigned* get_unsigned_for_string(const char* s) {
  unsigned v;
  int n = sscanf(s, "%u", &v);
  if(n != 1) return NULL;
  unsigned* pv = (unsigned*)malloc(sizeof(unsigned));
  if(pv == NULL) return NULL;
  *pv = v;
  return pv;
}


/**Sets a quikmat from a string.
 *\param m ptr to a quikmat struct
 *\param matstring string to convert from
 *\return non-zero if there was an error; zero if success
 */
int set_quikmat(struct quikmat* m, char* matstring) {
  char** matpieces;
  matpieces = g_strsplit(matstring, "/", 4);
  if(matpieces[0] == NULL) {
    /*This is a single material.*/
    GString* s = g_string_new(matstring);
    m->mat1 = s->str;
    g_string_free(s, FALSE);
    m->mat2 = NULL;
    m->x = 1.0;
    return 0;
  }
  GString* s = g_string_new(matpieces[0]);
  m->mat1 = s->str;
  g_string_free(s, FALSE);
  if(matpieces[1] == NULL) {
    /*Again a single material.  Ship it.*/
    m->mat2 = NULL;
    m->x = 1.0;
    return 0;
  }
  int val = sscanf(matpieces[1], "%lg", &(m->x));
  if((val != 1) || isnan(m->x) || isinf(m->x) || (m->x <=0) || (m->x >= 1)) {
    fprintf(stderr, "Got invalid number %g (%s) as material percentage in material string (%s)\n", m->x, matpieces[1], matstring);
    m->mat2 = NULL;
    return 1;
  }
  if(matpieces[2] == NULL) {
    fprintf(stderr, "Got no second alloy material in material string (%s)\n", matstring);
    m->mat2 = NULL;
    return 2;
  }
  m->mat2 = matpieces[2];
  return 0;
}

int process_args(int argc, char** argv, struct args* args) {
  int errored = 0;
  gchar **file_names = NULL;
  gboolean print_version = FALSE;
  gchar** parameters=NULL;
  GString *s = g_string_new("/home/solarion/dist/share/materials-databases/matdb.txt");
  args->matdb = s->str;
  args->N_starting_points=10;
  /*Set the random seed value to a random value (done internally when seed is 0) by default*/
  args->seed_value=0;
  args->blocksize=4;
  g_string_free(s, FALSE);
  s = g_string_new("OptimizerMaterial");
  args->optimizer_material = s->str;
  args->continuous = FALSE;
  #ifdef HAVE_MPI
  args->mpi = FALSE;
  args->temperature = 0;
  args->anneal_start = 0;
  args->anneal_end = 0;
  args->num_anneal_steps = 0;
  args->anneal_cycles = 0;
  args->calculation_steps_per_annealing_step = 0;
  #endif
  g_string_free(s, FALSE);
  char* mdb = NULL;
  char* optmat = NULL;
  int raw_startingpoints = 10;
  char* mat1=NULL;
  char* mat0=NULL;
  char* assignments = NULL;
  args->assignments = NULL;
  args->assignment_expr = NULL;
  int raw_blocksize;
  args->minx=0.0;
  args->maxx=1.0;
  args->stepsize=0.1;
  args->min_stepsize=0.001;
  args->multiplier=0.33;
  args->N_blocks = 0;
  args->do_nofield_comparison = FALSE;
  char* annealing_program=NULL;
  GOptionEntry options[] = {
    { G_OPTION_REMAINING, '\0', 0, G_OPTION_ARG_FILENAME_ARRAY, &file_names,},
    { "version", 'v', G_OPTION_FLAG_OPTIONAL_ARG, G_OPTION_ARG_NONE, &print_version, "Print version information about this program", NULL },
    { "matdb", 'm', G_OPTION_FLAG_OPTIONAL_ARG, G_OPTION_ARG_STRING, &(mdb), "Use an alternate materials database file.", NULL },
    {"parameter", 'p', G_OPTION_FLAG_OPTIONAL_ARG, G_OPTION_ARG_STRING_ARRAY, &parameters, "Provide a template parameter argument, in the form param,value (single value) or param,start,stop,step (sequence)\n", NULL},
    { "optmat", 'o', G_OPTION_FLAG_OPTIONAL_ARG, G_OPTION_ARG_STRING, &(optmat), "Material name specifying the region to be optimized", NULL },
    { "mat1", '1', G_OPTION_FLAG_OPTIONAL_ARG, G_OPTION_ARG_STRING, &(mat1), "Material name specifying one material to use in the optimizer.  Must be of the form mata/x/matb for alloys of MatA_x MatB_{1-x}.", NULL },
    { "mat0", '0', G_OPTION_FLAG_OPTIONAL_ARG, G_OPTION_ARG_STRING, &(mat0), "Material name specifying the other material to be used in the optimizer.  Must be of the form mata/x/matb for alloys of MatA_x MatB_{1-x}.", NULL },
    { "startingpoints", 's', G_OPTION_FLAG_OPTIONAL_ARG, G_OPTION_ARG_INT, &(raw_startingpoints), "Number of random starting structures", NULL },
    { "blocksize", 'b', G_OPTION_FLAG_OPTIONAL_ARG, G_OPTION_ARG_INT, &(raw_blocksize), "Size of a block of material in gridsites (optimization is done in blocks for accuracy; ignored if num_blocks is set to a nonzero value)", NULL },
    { "num_blocks", 'B', G_OPTION_FLAG_OPTIONAL_ARG, G_OPTION_ARG_INT, &(args->N_blocks), "Number of blocks to use (leave unspecified to use blocksize instead)", NULL },
    /*The seed can be signed here; this will allow full access to the bits*/
    { "seed", 'S', G_OPTION_FLAG_OPTIONAL_ARG, G_OPTION_ARG_INT64, (long unsigned*)&(args->seed_value), "Seed to use for the random number generator (to reproduce a run; default is random)", NULL },
    { "minx", 'x', G_OPTION_FLAG_OPTIONAL_ARG, G_OPTION_ARG_DOUBLE, &(args->minx), "Minimum allowable alloying for continuous alloying. (default is 0; ignored for digital alloys)", NULL},
    { "maxx", 'X', G_OPTION_FLAG_OPTIONAL_ARG, G_OPTION_ARG_DOUBLE, &(args->maxx), "Maximum allowable alloying for continuous alloying (default is 1; ignored for digital alloys)", NULL},
    { "multiplier", 'u', G_OPTION_FLAG_OPTIONAL_ARG, G_OPTION_ARG_DOUBLE, &(args->multiplier), "Alloying step multiplier (sets the next alloying stepsize after a successful optimization at one alloying stepsize) for continuous alloying (ignored for digital alloys)", NULL},
    { "stepsize", 't', G_OPTION_FLAG_OPTIONAL_ARG, G_OPTION_ARG_DOUBLE, &(args->stepsize), "Initial alloying step size (how much to change x for the first optimization step; subsequent steps are set by stepsize=multiplier*stepsize) for continuous alloying (ignored for digital alloys)", NULL},
    { "converged_stepsize", 'T', G_OPTION_FLAG_OPTIONAL_ARG, G_OPTION_ARG_DOUBLE, &(args->min_stepsize), "Consider alloying converged when stepsize reaches this value or smaller for continuous alloying (ignored for digital alloys)", NULL},
    {"assignments", 'a', G_OPTION_FLAG_OPTIONAL_ARG, G_OPTION_ARG_STRING, &assignments, "Provide initial material assignments.  Can be a single float, a comma-separated list of floats, or a function.  If digital alloying, 0 <= v < 0.5 is assigned 0; 0.5 <= v <= 1.0 is assigned 1.\n", NULL},
    { "continuous", 'c', G_OPTION_FLAG_OPTIONAL_ARG, G_OPTION_ARG_NONE, &(args->continuous), "Do a continuous optimization, not a digital alloying. NOTE that mat1 and mat2 must be a single material, not an alloy", NULL },
    #ifdef HAVE_MPI
    { "mpi", 'M', G_OPTION_FLAG_OPTIONAL_ARG, G_OPTION_ARG_NONE, &(args->mpi), "Use the MPI variant of the digital alloying algorithm.  Currently ignored for continuous optimization.  If an annealing option is not given, it uses a steepest-descent-related Monte Carlo method.", NULL },
    { "anneal-temperature", 'A', G_OPTION_FLAG_OPTIONAL_ARG, G_OPTION_ARG_DOUBLE, &(args->temperature), "If specified, simulated annealing will be used at the specified temperature instead of the steepest-descent mode.  Requires the MPI option.  Temperature is in units of 1e-13 eV m.", NULL },
    { "annealing-program", 'P', G_OPTION_FLAG_OPTIONAL_ARG, G_OPTION_ARG_STRING, &(annealing_program), "If specified, simulated annealing will be used with the specified program instead of the steepest-descent mode.  Requires the MPI option.  The program is specified in the form of num_cycles,start_temperature,stop_temperature,num_temperature_steps.  Temperature is in units of 1e-13 eV m.", NULL },
    #endif
    { "compare_noefield", '\0', G_OPTION_FLAG_OPTIONAL_ARG, G_OPTION_ARG_NONE, &(args->do_nofield_comparison), "Instead of simply optimizing alpha, optimize the *change* in alpha between the structure as specified and the structure with no applied electric field.", NULL },
    /*At some point, we should allow for a starting point*/
    {NULL}
  };
  GOptionContext *ctx;
  ctx = g_option_context_new("- sets up a structure for k-dot-p calculations");
  g_option_context_add_main_entries(ctx, options, "Options");

  g_option_context_parse(ctx, &argc, &argv, NULL);

  gchar* htext = g_option_context_get_help(ctx, FALSE, NULL);
  g_option_context_free(ctx);

  if(mdb != NULL) {
    free(args->matdb);
    args->matdb = mdb;
  }
  if(optmat != NULL) {
    free(args->optimizer_material);
    args->optimizer_material = optmat;
  }

  #ifdef HAVE_MPI
  if(args->temperature != 0) args->temperature *= 1e-13;
  if((args->temperature != 0) && (annealing_program != NULL)) {
    fprintf(stderr, "ERROR: both a fixed annealing temperature as well as an annealing program were specified!  Only one may be given at a time.\n");
    errored=1;
  }
  if(annealing_program != NULL) {
    char** segments;
    int N_temp_steps;
    segments = g_strsplit(annealing_program, ",", 5);
    if(segments == NULL) {
      fprintf(stderr, "ERROR: unable to break the specified annealing program (%s) into comma-delimited parts.\n", annealing_program);
      errored=1;
    }else if(segments[0] == NULL) {
      fprintf(stderr, "ERROR: unable to find the number of annealing cycles in the annealing program (%s).\n", annealing_program);
      errored=1;
    }else if(segments[1] == NULL) {
      fprintf(stderr, "ERROR: unable to find the starting temperature in the annealing program (%s).\n", annealing_program);
      errored=1;
    }else if(segments[2] == NULL) {
      fprintf(stderr, "ERROR: unable to find the stopping temperature in the annealing program (%s).\n", annealing_program);
      errored=1;
    }else if(segments[3] == NULL) {
      fprintf(stderr, "ERROR: unable to find the number of temperature steps per cycle in the annealing program (%s).\n", annealing_program);
      errored=1;
    }else if(segments[4] == NULL) {
      fprintf(stderr, "ERROR: unable to find the number of calculations per annealing step in the annealing program (%s).\n", annealing_program);
      errored=1;
    }else{
      /*All the values are there; calculate them*/
      double* dv;
      unsigned* iv;
      iv = get_unsigned_for_string(segments[0]);
      if(iv == NULL) {
  fprintf(stderr, "ERROR: unable to find an unsigned integer for the number of annealing cycles (%s).\n", segments[0]);
  errored=1;
      }else{
  args->anneal_cycles = *iv;
  if(*iv == 0) {
    fprintf(stderr, "ERROR: specifying 0 annealing cycles is useless; I'm assuming this was in error.\n");
    errored=1;
  }
  free(iv);
      }
      dv = get_double_for_string(segments[1]);
      if(dv == NULL) {
  fprintf(stderr, "ERROR: unable to find an double-precision number for the annealing cycle start (%s).\n", segments[1]);
  errored=1;
      }else{
  args->anneal_start = *dv*1e-13;
  if(args->anneal_start <= 0) {
    fprintf(stderr, "ERROR: annealing cycle starting temperature must be > 0 (got %g).\n", *dv);
    errored=1;
  }
  free(dv);
      }
      dv = get_double_for_string(segments[2]);
      if(dv == NULL) {
  fprintf(stderr, "ERROR: unable to find an double-precision number for the annealing cycle end (%s).\n", segments[2]);
  errored=1;
      }else{
  args->anneal_end = *dv*1e-13;
  if(args->anneal_start <= 0) {
    fprintf(stderr, "ERROR: annealing cycle ending temperature must be > 0 (got %g).\n", *dv);
    errored=1;
  }
  if(args->anneal_start >= args->anneal_end) {
    fprintf(stderr, "ERROR: annealing starting temperature (%g) is greater than or equal to the ending temperature (%g)!\n", args->anneal_start, args->anneal_end);
    errored=1;
  }
  free(dv);
      }
      iv = get_unsigned_for_string(segments[3]);
      if(iv == NULL) {
  fprintf(stderr, "ERROR: unable to find an unsigned integer for the number of annealing steps per cycle (%s).\n", segments[3]);
  errored=1;
      }else{
  args->num_anneal_steps = *iv;
  if(*iv == 0) {
    fprintf(stderr, "ERROR: specifying 0 annealing steps is useless; I'm assuming this was in error.\n");
    errored=1;
  }
  free(iv);
      }
      iv = get_unsigned_for_string(segments[4]);
      if(iv == NULL) {
  fprintf(stderr, "ERROR: unable to find an unsigned integer for the number of calculation steps per annealing step (%s).\n", segments[4]);
  errored=1;
      }else{
  args->calculation_steps_per_annealing_step = *iv;
  if(*iv == 0) {
    fprintf(stderr, "ERROR: specifying 0 calculations per annealing step is useless; I'm assuming this was in error.\n");
    errored=1;
  }
  free(iv);
      }
    }
  }
  #endif

  if(raw_startingpoints <= 0) {
    fprintf(stderr, "ERROR: starting points given (%d) was less than or equal to zero! This is nonsensical.\n", raw_startingpoints);
    errored=1;
  }else{
    args->N_starting_points=raw_startingpoints;
  }
  if(raw_blocksize <= 0) {
    fprintf(stderr, "ERROR: block size (%d) given was less than or equal to zero! This is nonsensical.\n", raw_blocksize);
    errored=1;
  }else{
    args->blocksize=raw_blocksize;
  }

  if((mat1 == NULL) || (mat0 == NULL)) {
    fprintf(stderr, "ERROR: you MUST supply the two material types for the digital alloy optimizer.\n");
    errored=1;
  }else{
    fprintf(stderr, "mat1=(%s) mat0=(%s)\n", mat1, mat0);
    errored=set_quikmat(&(args->mat1), mat1);
    errored=set_quikmat(&(args->mat0), mat0);
  }

  args->filestreams = NULL;
  GFileInputStream *instream;
  GError *gerr = NULL;
  if(file_names != NULL) {
    for(int i = 0; file_names[i] != NULL; i++) {
      struct file_info *fi = (struct file_info*)malloc(sizeof(struct file_info));
      if(fi == NULL) {
  fprintf(stderr, "ERROR: unable to allocate input file info struct!\n");
  return 500;
      }
      if((file_names[i][0] == '-') && (file_names[i][1] == '\0')) {
  fi->instream = (GInputStream*)g_unix_input_stream_new(0, FALSE);
  fi->filename = g_string_new("-");
      }else{
  GFile *f = g_file_new_for_commandline_arg(file_names[i]);
  gerr = NULL;
  instream = g_file_read(f, NULL, &gerr);
  if(gerr != NULL) {
    fprintf(stderr, "makestruct: ERROR reading in structure file (%s): %s", g_file_get_uri(f), gerr->message);
    return 100;
  }
  g_object_unref(f);
  fi->instream = (GInputStream*)instream;
  fi->filename = g_string_new(file_names[i]);
      }
      args->filestreams = g_list_append(args->filestreams, fi);
    }
  }
  if(args->filestreams == NULL) {
    fprintf(stderr, "ERROR: you MUST give at least one file to process!\n");
    errored=1;
  }

  if(args->continuous) {
    if(args->mat1.x != 1) {
      fprintf(stderr, "ERROR: you MUST not give an alloy for continuous alloying (material 1 was an alloy.)\n");
      errored=1;
    }
    if(args->mat0.x != 1) {
      fprintf(stderr, "ERROR: you MUST not give an alloy for continuous alloying (material 0 was an alloy.)\n");
      errored=1;
    }
  }
  
  if(assignments != NULL) {
    char** strs = g_strsplit(assignments, ",", -1);
    if(strs[0] == NULL) {
      fprintf(stderr, "ERROR: initial material assignment argument was given, but no data was!\n");
      errored == 1;
    }else{
      /*Try to convert the first argument into a double.  If we can't, it's an equation.*/
      double* pv = get_double_for_string(strs[0]);
      if(pv == NULL) {
  if(strs[1] != NULL) {
    /*Need to set up the equation parser*/
    fprintf(stderr, "ERROR: First number given was not a float (%s), but you provided multiple arguments!  You may only specify one equation to use for initial material assignments.\n", strs[0]);
    errored == 1;
  }else{
    /*Need to set up the equation parser*/
    fprintf(stderr, "ERROR: Equations are not yet supported.  Sorry.  One thing at a time.\n");
    errored == 1;
  }
      }
       /*If the first arg wasn't an equation, we're good.*/
      if(pv != NULL) {
  free(pv);
  for(unsigned i=0; strs[i] != NULL; i++) {
    double* pv = get_double_for_string(strs[i]);
    if(pv == NULL) {
      fprintf(stderr, "ERROR: Unable to convert %udth number (%s) to a float!\n", i, strs[i]);
      errored == 1;
      break;
    }
    /*OK.  I'm a softie.  <0.5 => 0; >=0.5 => 1*/
    if(args->continuous) {
      if(*pv < 0.5) {
        *pv = 0.0;
      }else{
        *pv = 1.0;
      }
    }
    args->assignments = g_list_append(args->assignments, pv);
  }
      }
    }
    g_strfreev(strs);
  }
  
  char** tparam_entry;
  args->template_parameters = NULL;
  if(parameters != NULL) {
    for(unsigned i=0; parameters[i] != NULL; i++) {
      tparam_entry = g_strsplit(parameters[i], ",", -1);
      int len=g_strv_length(tparam_entry);
      struct template_parameter_range* tpr;
      if(len == 2) {
  /*parameter,value form*/
  tpr = new_template_parameter_range(tparam_entry[0], strtod(tparam_entry[1], NULL), strtod(tparam_entry[1], NULL), strtod(tparam_entry[1], NULL));
      }else if(len == 4) {
  tpr = new_template_parameter_range(tparam_entry[0], strtod(tparam_entry[1], NULL), strtod(tparam_entry[2], NULL), strtod(tparam_entry[3], NULL));
      }else{
  fprintf(stderr, "ERROR: invalid parameteter entry: incorrect number of parameter range arguments (got %d; expected 2 or 4; should be either parameter,value or parameter,start,stop,step)", len);
  return 15;
      }
      if(tpr == NULL) {
  fprintf(stderr, "ERROR: unable to allocate space for new template parameter range (parameter is %s)\n", tparam_entry[0]);
  return 16;
      }
      args->template_parameters = g_list_append(args->template_parameters, tpr);
    }
  }

  if(print_version){
    print_version_info();
    return 1;
  }
  if(errored){
    printf("%s", htext);
    printf("\n\n**For version and compilation information, run with --version***\n");
    return 3;
  }
  return 0;
}