//package:  starlab
//date:     11-feb-95
//author:   Piet Hut
//task:     mkcube
//does:     construct a simple homogeneous cube
//index:    model, init
//also:     mkplummer
//+desc
//	      M = 1, T/U = -1/2, E = -1/4
//  If the "-u" flag is set, the particles are left unscaled.
//-desc

#include "dyn.h"

local void  mkcube(dyn * root, int n, int u_flag)
{
    dyn  * bi;

    root->set_mass(1);
    real pmass = 1.0 / n;

    for (bi = root->get_oldest_daughter(); bi != NULL;
         bi = bi->get_younger_sister()) {
	bi->set_mass(pmass);	
        bi->set_pos(vector(randinter(-1,1),randinter(-1,1),randinter(-1,1)));
        bi->set_vel(vector(randinter(-1,1),randinter(-1,1),randinter(-1,1)));
    }

//  Transform to center-of-mass coordinates, and optionally
//  scale to standard parameters.

    root->to_com();

    if (! u_flag && n > 1) {

        real kinetic, potential;

	get_energies(root, 0.0, kinetic, potential);
	scale_virial(root, -0.5, kinetic, potential);
	scale_energy(root, -0.25, kinetic + potential);
    }
}

#ifdef TOOLBOX

#define  SEED_STRING_LENGTH  60

#  define  FALSE  0
#  define  TRUE   1

#define USAGE "usage: mkcube -n # [-s #] [-r #] [-m #] [-i] [-u] [-c \"...\"]\n"

//usage: mkcube -n # [-s #] [-r #] [-m #] [-i] [-u] [-c "..."]

main(int argc, char ** argv) {
    int  i;
    int  c;
    int  n;
    int  input_seed, actual_seed;
    int  c_flag = FALSE;
    int  n_flag = FALSE;
    int  s_flag = FALSE;
    int  i_flag = FALSE;
    int  l_flag = FALSE;
    int  u_flag = FALSE;   // if TRUE, leave scaling Untouched
    char  *comment;
    char  seedlog[SEED_STRING_LENGTH];
    extern char *poptarg;
    int  pgetopt(int, char **, char *);
    real  cube_size = 2;

    if (argc == 1) {
	cerr << USAGE;
	exit(1);
    }

    while ((c = pgetopt(argc, argv, "c:n:s:m:r:uil")) != -1)
	switch(c)
	    {
//option:   -n #
//          number of bodies

	    case 'n': n_flag = TRUE;
		      n = atoi(poptarg);
		      break;
//option:   -s #
//          Initial random number seed
	    case 's': s_flag = TRUE;
		      input_seed = atoi(poptarg);
		      break;
//option:   -i
//          Adds numeric indices to particles
	    case 'i': i_flag = TRUE;
		      break;
//option:   -l
	    case 'l': l_flag = TRUE;
		      break;
//option:   -u
	    case 'u': u_flag = TRUE;
		      break;
//option:   -c "..."
//          comment
	    case 'c': c_flag = TRUE;
		      comment = poptarg;
		      break;
            case '?': cerr << USAGE;
		      exit(1);
	    }            
    
    if (n_flag == FALSE) {
        cerr << "makecube: must specify the number # of";
	cerr << " particles with -n#\n";
	exit(1);
    }
    
    if (n < 1) {
        cerr << "mkcube: n < 1 not allowed\n";
	exit(1);
    }

    dyn *b, *by, *bo;
    b = new dyn();
    bo = new dyn();
    if (i_flag) 
        bo->set_label(1);
    b->set_oldest_daughter(bo);
    bo->set_parent(b);
    if (l_flag)
	putrq(b->get_dyn_story(), "cube_size", cube_size);

    for (i = 1; i < n; i++) {
        by = new dyn();
	if (i_flag) by->set_label(i+1);
	by->set_parent(b);
        bo->set_younger_sister(by);
        by->set_elder_sister(bo);
        bo = by;
    }

    if (c_flag == TRUE) b->log_comment(comment);

    b->log_history(argc, argv);

    if (s_flag == FALSE) input_seed = 0;
    actual_seed = srandinter(input_seed);

    sprintf(seedlog, "         random number generator seed = %d",actual_seed);
    b->log_comment(seedlog);

    if (n > 0) mkcube(b, n, u_flag);

    put_node(cout, *b);
}

#endif

/* end of: mkcube.c */