/*****************************************************************************
 * Copyright 2005 Kevin Hobbs <hobbsk@ohiou.edu>                             *
 *                                                                           *
 * This file is part of Balloon.                                             *
 *                                                                           *
 * Balloon is free software; you can redistribute it and/or modify           *
 * it under the terms of the GNU General Public License as published by      *
 * the Free Software Foundation; either version 2 of the License, or         *
 * (at your option) any later version.                                       *
 *                                                                           *
 * Balloon 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 General Public License for more details.                              *
 *                                                                           *
 * You should have received a copy of the GNU General Public License         *
 * along with Balloon; if not, write to the Free Software                    *
 * Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA *
 *****************************************************************************/


#include <stdio.h>
#include <stddef.h>
#include <stdlib.h>
#include "mesh.h"
#include "utility.h"
#include "balloon.h"


struct point ** 
init_balloon(
		const double * x0, 
		const double * s, 
		const size_t ndims,
		size_t * id,
	        gsl_rng * r	
	    )
{
	size_t i, j;
	struct point ** balloon;
	double * si;
 	
	balloon = xmalloc( sizeof( struct point * ) * ( ndims * 2 ) );
	
	/***********************************/
	/* Create N * 2 mesh points at x_0 */
	/***********************************/
	for ( i = 0; i < ndims * 2; i++)
	{
		balloon[i] = new_empty_point(ndims, *id);
		*id = (*id) + 1;
		for ( j = 0; j < ndims; j++)
			balloon[i] -> x_in[j] = x0[j];
	}
	
	/*************************************************/
	/* Each point should be connected to every other */
	/* except for the one on the same axis           */
	/*************************************************/
	for ( i = 0; i < ndims; i++)
		for ( j = 0; j < ndims; j++)
			if (i != j)
			{
				connect_points(balloon[i        ], balloon[j        ]);
				connect_points(balloon[i + ndims], balloon[j + ndims]);
				connect_points(balloon[i        ], balloon[j + ndims]);
				connect_points(balloon[i + ndims], balloon[j        ]);
			}

	/*************************************************/
	/* Push the leg points away along each dimension */
	/*************************************************/
	si = xmalloc( sizeof( double ) * ndims );
	for ( i = 0; i < ndims; i++)
	{
		for ( j = 0; j < ndims; j++)
			si[j] = 0;
		si[i] = s[i];
		inflate(balloon[i], si, s, ndims, r);
		si[i] = -s[i];
		inflate(balloon[i + ndims], si, s, ndims, r);
	}
	
		
	free(si);
	return balloon;
}