#include "peg_pic.h"


/*******  Allocate space for main parameters  ********/
/* 
 * WARNING! Anything that is allocated here, must have counterpart in freespace() - to deallocate  (free) memory 
 */

void alsp() {
	int i, isp, k, l;

	x_array= (float *)malloc(ng*sizeof(float));

	for(i=0; i< ng; i++) {
		x_array[i]= i*dx;
	}

#ifdef REFLECTED_OWN_SPECIES
	/* we have twice as much species because of particle reflection */
	nsp *= 2;
#endif

	/** Allocate space for general parameters **/

	phi= (float *)malloc(ng*sizeof(float));
	rho= (float *)malloc(ng*sizeof(float));
	e  = (float *)malloc(ng*sizeof(float));

	fr_E  = (float *)malloc(nsp*sizeof(float));
	jwall = (float *)malloc(nsp*sizeof(float));
	qm = (float *)malloc(nsp*sizeof(float));
	q  = (float *)malloc(nsp*sizeof(float));
	m  = (float *)malloc(nsp*sizeof(float));
	Escale = (float *)malloc(nsp*sizeof(float));
	Tscale = (float *)malloc(nsp*sizeof(float));
	Fscale = (float *)malloc(nsp*sizeof(float));
	np   = (int *)malloc(nsp*sizeof(int));
	nstep= (int *)malloc(nsp*sizeof(int));
	wt= (int *)malloc(nsp*sizeof(int));
	ht_sp= (int *)malloc(nsp*sizeof(int));
	dinj = (int *)malloc(nsp*sizeof(int));
	i_t  = (int *)malloc(nsp*sizeof(int));
	inject = (int *)malloc(nsp*sizeof(int));
	F_temp = (float *)malloc(nsp*sizeof(float));
	F_tempRHS = (float *)malloc(nsp*sizeof(float));
	FLHS = (float *)malloc(nsp*sizeof(float));
	FRHS = (float *)malloc(nsp*sizeof(float));

	gridref_r = (float *)malloc(nsp*sizeof(float));			/* current reflection on RIGHT grid border */
	gridref_r0 = (float *)malloc(nsp*sizeof(float));		/* initial reflection on RIGHT grid border */
	gridref_rsramp = (float *)malloc(nsp*sizeof(float));	/* size step of reflection on RIGHT grid border - values is multiplied by time */
	gridref_l = (float *)malloc(nsp*sizeof(float));			/* current reflection on LEFT grid border */
	gridref_l0 = (float *)malloc(nsp*sizeof(float));		/* initial reflection on LEFT grid border */
	gridref_lsramp = (float *)malloc(nsp*sizeof(float));	/* size step of reflection on LEFT grid border - values is multiplied by time */

	for(i=0; i < nsp; i++) {
		/* initial fill of grid reflection on left and right boundary */
		gridref_r[i] = 0.0;
		gridref_r0[i] = 0.0;
		gridref_rsramp[i] = 0.0;
		gridref_l[i] = 0.0;
		gridref_l0[i] = 0.0;
		gridref_lsramp[i] = 0.0;
	}

	Nsec = (int *)malloc(nsp*sizeof(int));
	srho = (float **)malloc(nsp*sizeof(float *));

	for(isp=0; isp<nsp; isp++) {
		srho[isp] = (float *)malloc(ng*sizeof(float));
	}


	/* Injection from walls */
	enter  = (float **)malloc(nsp*sizeof(float *));
	vnpart = (float **)malloc(nsp*sizeof(float *));
	vnperpt= (float **)malloc(nsp*sizeof(float *));
	vnpar0 = (float **)malloc(nsp*sizeof(float *));
	vnperp0= (float **)malloc(nsp*sizeof(float *));
	parvell= (float **)malloc(nsp*sizeof(float *));
	parvelr= (float **)malloc(nsp*sizeof(float *));

	for (isp=0; isp < nsp; isp++) 	{
		enter[isp]  = (float *)malloc(2*sizeof(float));
		vnpart[isp] = (float *)malloc(2*sizeof(float));
		vnperpt[isp]= (float *)malloc(2*sizeof(float));
		vnpar0[isp] = (float *)malloc(2*sizeof(float));
		vnperp0[isp]= (float *)malloc(2*sizeof(float));
		parvell[isp]= (float *)malloc(VVMAX*sizeof(float));
		parvelr[isp]= (float *)malloc(VVMAX*sizeof(float));
	}

	/* Initial loading */
	vpart  = (float *)malloc(nsp*sizeof(float));
	vpar0  = (float *)malloc(nsp*sizeof(float));
	vperpt = (float *)malloc(nsp*sizeof(float));
	vperp0 = (float *)malloc(nsp*sizeof(float));


	if (sourcD) {                /* Particle source */
		vspart = (float *)malloc(nsp*sizeof(float));
		vspar0 = (float *)malloc(nsp*sizeof(float));
		vsperpt= (float *)malloc(nsp*sizeof(float));
		vsperp0= (float *)malloc(nsp*sizeof(float));
		Src_con = (float *)malloc(nsp*sizeof(float));
		Stemp   = (float *)malloc(nsp*sizeof(float));
		
		if (ht>0.) {               /* "Heating" of the source */
			heat = (float *)malloc(nsp*sizeof(float));
		}
		
		if( sourcD == 4 ){
			/* allocate additional arrays... */
			totalLen = (double *)malloc(ng*sizeof(double)) ;
			kLen = (double *)malloc(ng*sizeof(double)) ;
			phi_six = (struct pt *)malloc(ng*sizeof(struct pt));
			
			for(i=0; i< ng; i++) {
				phi_six[i].x = i;
				phi_six[i].y = 0.0;
			}
		}
		
	} else if (ccc) {  /* For calculation of Landau log. */
		vspart = (float *)malloc(nsp*sizeof(float));
	}
	
	if (sourcD1) {           /* Second particle source */
		vspart1 = (float *)malloc(nsp*sizeof(float));
		vspar01 = (float *)malloc(nsp*sizeof(float));
		vsperpt1= (float *)malloc(nsp*sizeof(float));
		vsperp01= (float *)malloc(nsp*sizeof(float));
		Src_con1 = (float *)malloc(nsp*sizeof(float));
		Stemp1   = (float *)malloc(nsp*sizeof(float));

		if (ht1>0.) {              /* "Heating" of the source */
			heat1 = (float *)malloc(nsp*sizeof(float));
		}
	}


	if (field) {
		av_e = (float **)malloc(nsp*sizeof(float *));
		
		for (isp=0; isp < nsp; isp++) {
			av_e[isp] = (float *)malloc(ng*sizeof(float));
		}
	}


	if (Ey_flag) {               /* Field parallel to the wall */
		Ey = (float *)malloc(nsp*sizeof(float));
		E1y= (float *)malloc(nsp*sizeof(float));
	}

	if (bmag > 0.0)	{
		tkv= (float *)malloc(nsp*sizeof(float));
		tx = (float *)malloc(nsp*sizeof(float));
		tz = (float *)malloc(nsp*sizeof(float));
		sx = (float *)malloc(nsp*sizeof(float));
		sz = (float *)malloc(nsp*sizeof(float));
		count_max = (int *)malloc(nsp*sizeof(int));
	}


	if (ccc) {              /* Coulomb collisions */
		Ldav = (int ***)malloc(ccc*sizeof(int **));
		pp = (float *)malloc(ccc*sizeof(float));
		LL = (float *)malloc(ccc*sizeof(float));
		Z  = (float *)malloc(ccc*sizeof(float));
		miu= (float *)malloc(ccc*sizeof(float));

		for(isp=0; isp<=ccc; isp++) {
			Ldav[isp] = (int **)malloc(nc*sizeof(int *));
		}

		for(isp=0; isp<=ccc; isp++) {
			for (i=0; i< nc; i++) {
				Ldav[isp][i] = (int *)malloc(MAXD*sizeof(int));
			}
		}
	}

	if (Nnum) {           /* Charged-neutral particle collisions */
		Nneut = (float **)malloc(Nnum*sizeof(float *));
		Tneut = (float **)malloc(Nnum*sizeof(float *));

		if (slow) {
			Tneutshow = (float **)malloc(Nnum*sizeof(float *));
		}

		for (i=0; i < Nnum; i++) {
			Nneut[i] = (float *)malloc(ng*sizeof(float));
			Tneut[i] = (float *)malloc(ng*sizeof(float));
			if (slow) {
				Tneutshow[i] = (float *)malloc(ng*sizeof(float));
			}
		}

		if ((ecollisional) && (Nnum)) {
			sctane= (float **)malloc(3*sizeof(float *)); // For scat. angle
			Nexel = (int *)malloc(Nnum*sizeof(int)); // Max. number of excit. coll.

			for (i=0; i<3; i++) {
				sctane[i]= (float *)malloc(Nnum*sizeof(float));
			}
		}

		k= 1;
		
		for (isp=0; isp<Nnum; isp++) {
			for (i=0; i < Nnum; i++) {
				if (icollisional[isp][i] && k) {
					nimiu = (float **)malloc(Nnum*sizeof(float *));
					sctang= (float **)malloc(Nnum*sizeof(float *));
					
					for (l=0; l < Nnum; l++) {
						nimiu[l] = (float *)malloc(Nnum*sizeof(float));
						sctang[l]= (float *)malloc(Nnum*sizeof(float));
					}
					
					k= 0;
				}
			}
		}
	}  /* end of charged-neutral particle collisions */

	/* Surface physics */

	srf = (int **)malloc(nsp*sizeof(int *));
	secsp = (int **)malloc(nsp*sizeof(int *));
	sectype = (int **)malloc(nsp*sizeof(int *));

	injtype= (int ***)malloc(nsp*sizeof(int **));
	seec  = (float ***)malloc(nsp*sizeof(float **));
	Tseec = (float ***)malloc(nsp*sizeof(float **));
	Ew    = (float ***)malloc(nsp*sizeof(float **));
	E0    = (float ***)malloc(nsp*sizeof(float **));

	for(isp=0; isp<nsp; isp++) {
		srf[isp]= (int *)malloc(2*sizeof(int));
		secsp[isp]= (int *)malloc(MAXSRF*sizeof(int));
		sectype[isp]= (int *)malloc(MAXSRF*sizeof(int));

		injtype[isp]= (int **)malloc(MAXSRF*sizeof(int *));
		seec[isp] = (float **)malloc(MAXSRF*sizeof(float *));
		Tseec[isp]= (float **)malloc(MAXSRF*sizeof(float *));
		Ew[isp]   = (float **)malloc(MAXSRF*sizeof(float *));
		E0[isp]   = (float **)malloc(MAXSRF*sizeof(float *));
	}

	for(isp=0; isp<nsp; isp++) {
		for(i=0; i<MAXSRF; i++) {
			injtype[isp][i]= (int *)malloc(2*sizeof(int));
			seec[isp][i] = (float *)malloc(2*sizeof(float));
			Tseec[isp][i]= (float *)malloc(2*sizeof(float));
			Ew[isp][i]   = (float *)malloc(2*sizeof(float));
			E0[isp][i]   = (float *)malloc(2*sizeof(float));
		}
	}
	
	if (slow || stopfl) {
		/** Allocate space for diagnostics **/
		Vx= (float **)malloc(nsp*sizeof(float *));
		Vpar= (float **)malloc(nsp*sizeof(float *));
		Vy= (float **)malloc(nsp*sizeof(float *));
		Vz= (float **)malloc(nsp*sizeof(float *));
		Tx= (float **)malloc(nsp*sizeof(float *));
		Qh= (float **)malloc(nsp*sizeof(float *));
		Qtot= (float **)malloc(nsp*sizeof(float *));
		T_dif= (float **)malloc(nsp*sizeof(float *));
		Qscale= (float *)malloc(nsp*sizeof(float));
		Qwscale = (float *)malloc(nsp*sizeof(float));

		nbin = (int *)malloc(nsp*sizeof(int));
		nbin_fvx = (int *)malloc(nsp*sizeof(int));

		emin   = (float *)malloc(nsp*sizeof(float));
		de     = (float *)malloc(nsp*sizeof(float));
		dtheta = (float *)malloc(nsp*sizeof(float));
		Eflhs  = (float *)malloc(nsp*sizeof(float));
		Eflhs_temp= (float *)malloc(nsp*sizeof(float));
		Efrhs  = (float *)malloc(nsp*sizeof(float));
		Efrhs_temp = (float *)malloc(nsp*sizeof(float));
		Tlhs   = (float *)malloc(nsp*sizeof(float));
		Tlhs_temp= (float *)malloc(nsp*sizeof(float));
		Trhs   = (float *)malloc(nsp*sizeof(float));
		Trhs_temp = (float *)malloc(nsp*sizeof(float));
		de_fvx = (float *)malloc(nsp*sizeof(float));
		vx_min = (float *)malloc(nsp*sizeof(float));

		for(isp=0; isp<nsp; isp++) {
			Vx[isp] = (float *)malloc(ng*sizeof(float));
			Vpar[isp] = (float *)malloc(ng*sizeof(float));
			Vy[isp] = (float *)malloc(ng*sizeof(float));
			Vz[isp] = (float *)malloc(ng*sizeof(float));
			Tx[isp] = (float *)malloc(ng*sizeof(float));
			Qh[isp] = (float *)malloc(ng*sizeof(float));
			Qtot[isp] = (float *)malloc(ng*sizeof(float));
			T_dif[isp] = (float *)malloc(ng*sizeof(float));
		}

		/** Allocate space for collision diagnostics **/

		if (ecollisional) {
			ionrate  = (float **)malloc(Nnum*sizeof(float *));
			extrate  = (float **)malloc(Nnum*sizeof(float *));

			for (isp=0; isp< Nnum; isp++) {
				ionrate[isp]  = (float *)malloc(ng*sizeof(float));
				extrate[isp]  = (float *)malloc(ng*sizeof(float));
				
				for (i=0; i<ng; i++) {
					ionrate[isp][i]= extrate[isp][i]= 0;
				}
			}
		}

		if (icollns) {
			chrgxrate= (float **)malloc(Nnum*sizeof(float *));
			
			for (isp=0; isp < Nnum; isp++) {
				chrgxrate[isp]= (float *)malloc(ng*sizeof(float));
				
				for (i=0; i<ng; i++) {
					chrgxrate[isp][i]= 0;
				}
			}
		}

	/* end of (slow || stopfl) */
	} else {
		if (Ef_flag) {        /* for energy flux */
			Eflhs  = (float *)malloc(nsp*sizeof(float));
			Eflhs_temp= (float *)malloc(nsp*sizeof(float));
			Efrhs  = (float *)malloc(nsp*sizeof(float));
			Efrhs_temp = (float *)malloc(nsp*sizeof(float));
			Tlhs   = (float *)malloc(nsp*sizeof(float));
			Tlhs_temp= (float *)malloc(nsp*sizeof(float));
			Trhs   = (float *)malloc(nsp*sizeof(float));
			Trhs_temp = (float *)malloc(nsp*sizeof(float));
			Qwscale = (float *)malloc(nsp*sizeof(float));
		}

		if (qtp>2) {          /* for effective ionization and/or movie */
			Vxef= (float *)malloc(ng*sizeof(float));
			Vyef= (float *)malloc(ng*sizeof(float));
			Vzef= (float *)malloc(ng*sizeof(float));
			Tx= (float **)malloc(nsp*sizeof(float *));

			if (qtp==3) {
				Vpar= (float **)malloc(nsp*sizeof(float *));
				Qtot= (float **)malloc(nsp*sizeof(float *));
				Qh= (float **)malloc(nsp*sizeof(float *));
				Qscale= (float *)malloc(nsp*sizeof(float));
			}

			for(isp=0; isp<nsp; isp++) {
				Tx[isp] = (float *)malloc(ng*sizeof(float));
				
				if (qtp==3)	{
					Vpar[isp] = (float *)malloc(ng*sizeof(float));
					Qtot[isp] = (float *)malloc(ng*sizeof(float));
					Qh[isp] = (float *)malloc(ng*sizeof(float));
				}
			}
		}

		if (!efi)  efiFL= 0;      /* no effective ionization */

	}   /* end of else */


	/** Allocate space for time history diagnostics **/
	t_array= (float *)malloc(HISTMAX*sizeof(float));
	np_hist= (float **)malloc(nsp*sizeof(float *));
	com_cur_hist= (float *)malloc(HISTMAX*sizeof(float));
	com_phi_hist= (float **)malloc(2*sizeof(float *));
	com_phi_hist[0]= (float *)malloc(HISTMAX*sizeof(float));
	com_phi_hist[1]= (float *)malloc(HISTMAX*sizeof(float));
	F_hist= (float **)malloc(nsp*sizeof(float *));
	F_histRHS= (float **)malloc(nsp*sizeof(float *));
	gridref_r_hist = (float **)malloc(nsp*sizeof(float));
	gridref_l_hist = (float **)malloc(nsp*sizeof(float));

	if (slow || stopfl || Ef_flag) {
		Eflhs_hist=(float **)malloc(nsp*sizeof(float *));
		Efrhs_hist=(float **)malloc(nsp*sizeof(float *));
		Tlhs_hist =(float **)malloc(nsp*sizeof(float *));
		Trhs_hist =(float **)malloc(nsp*sizeof(float *));
		Halhs =(float *)malloc(Nnum*sizeof(float));
		Harhs =(float *)malloc(Nnum*sizeof(float));
		Halhs_temp=(float *)malloc(Nnum*sizeof(float));
		Harhs_temp=(float *)malloc(Nnum*sizeof(float));
		Halhs_hist=(float **)malloc(Nnum*sizeof(float *));
		Harhs_hist=(float **)malloc(Nnum*sizeof(float *));

		for (i=0; i < Nnum; i++) {
			Halhs_hist[i] =(float *)malloc(HISTMAX*sizeof(float));
			Harhs_hist[i] =(float *)malloc(HISTMAX*sizeof(float));
		}
	}

	for(i=0; i<nsp; i++) {
		np_hist[i]= (float *)malloc(HISTMAX*sizeof(float));
		F_hist[i]= (float *)malloc(HISTMAX*sizeof(float));
		F_histRHS[i]= (float *)malloc(HISTMAX*sizeof(float)); 
		gridref_r_hist[i] = (float *)malloc(HISTMAX*sizeof(float));
		gridref_l_hist[i] = (float *)malloc(HISTMAX*sizeof(float));

		if (slow || stopfl || Ef_flag) {
			Eflhs_hist[i]= (float *)malloc(HISTMAX*sizeof(float));
			Efrhs_hist[i]= (float *)malloc(HISTMAX*sizeof(float));
			Tlhs_hist[i]= (float *)malloc(HISTMAX*sizeof(float));
			Trhs_hist[i]= (float *)malloc(HISTMAX*sizeof(float));
		}
	}

	if (slow || stopfl) {
		ese_hist= (float *)malloc(HISTMAX*sizeof(float));
		Teq = (float *)malloc(HISTMAX*sizeof(float));
		kes_hist= (float **)malloc(nsp*sizeof(float *));
		jwall_hist= (float **)malloc(nsp*sizeof(float *));
		Tx_hist=(float **)malloc(nsp*sizeof(float *));
		T_dif_hist=(float **)malloc(nsp*sizeof(float *));
		Vxt=(float **)malloc(nsp*sizeof(float *));

		for(i=0; i<nsp; i++) {
			kes_hist[i]= (float *)malloc(HISTMAX*sizeof(float));
			jwall_hist[i]= (float *)malloc(HISTMAX*sizeof(float));
			Tx_hist[i]= (float *)malloc(HISTMAX*sizeof(float));
			T_dif_hist[i]= (float *)malloc(HISTMAX*sizeof(float));
			Vxt[i]= (float *)malloc(HISTMAX*sizeof(float));
		}
	}

	/** Allocate space for frequency history diagnostics **/
	if(nfft) {
		srho_ave= (float **)malloc(nsp*sizeof(float *));
		srho_ave_show= (float **)malloc(nsp*sizeof(float *));
		Vx_ave= (float **)malloc(nsp*sizeof(float *));
		Vx_ave_show= (float **)malloc(nsp*sizeof(float *));
		Vpar_ave= (float **)malloc(nsp*sizeof(float *));
		Vpar_ave_show= (float **)malloc(nsp*sizeof(float *));
		Tx_ave= (float **)malloc(nsp*sizeof(float *));
		Qh_ave= (float **)malloc(nsp*sizeof(float *));
		Qtot_ave= (float **)malloc(nsp*sizeof(float *));
		Tx_ave_show= (float **)malloc(nsp*sizeof(float *));
		Qh_ave_show= (float **)malloc(nsp*sizeof(float *));
		Qtot_ave_show= (float **)malloc(nsp*sizeof(float *));
		T_dif_ave= (float **)malloc(nsp*sizeof(float *));
		T_dif_ave_show= (float **)malloc(nsp*sizeof(float *));

		for(isp=0; isp<nsp; isp++) {
			srho_ave[isp] = (float *)malloc(ng*sizeof(float));
			srho_ave_show[isp] = (float *)malloc(ng*sizeof(float));
			Vx_ave[isp] = (float *)malloc(ng*sizeof(float));
			Vx_ave_show[isp] = (float *)malloc(ng*sizeof(float));
			Vpar_ave[isp] = (float *)malloc(ng*sizeof(float));
			Vpar_ave_show[isp] = (float *)malloc(ng*sizeof(float));
			Tx_ave[isp] = (float *)malloc(ng*sizeof(float));
			Qh_ave[isp] = (float *)malloc(ng*sizeof(float));
			Qtot_ave[isp] = (float *)malloc(ng*sizeof(float));
			Tx_ave_show[isp] = (float *)malloc(ng*sizeof(float));
			Qh_ave_show[isp] = (float *)malloc(ng*sizeof(float));
			Qtot_ave_show[isp] = (float *)malloc(ng*sizeof(float));
			T_dif_ave[isp] = (float *)malloc(ng*sizeof(float));
			T_dif_ave_show[isp] = (float *)malloc(ng*sizeof(float));

			for(i=0; i< ng; i++) {
				srho_ave[isp][i]= srho_ave_show[isp][i]= Vx_ave[isp][i]= 0.0;
				Vx_ave_show[isp][i] = Tx_ave[isp][i] = Tx_ave_show[isp][i] =0.0;
				T_dif_ave[isp][i] =T_dif_ave_show[isp][i] = Vpar[isp][i]= 0.0;
				Vpar_ave[isp][i]=Vpar_ave_show[isp][i]= Qh_ave[isp][i] = 0.0;
				Qh_ave_show[isp][i]= Qtot_ave[isp][i]= Qtot_ave_show[isp][i]= 0.0;
			}
		}

		phi_ave= (float *)malloc(ng*sizeof(float));
		phi_ave_show= (float *)malloc(ng*sizeof(float));
		
		for(i=0; i< ng; i++) {
			phi_ave[i]= phi_ave_show[i]= 0.0;
		}

		e_ave= (float *)malloc(ng*sizeof(float));
		e_ave_show= (float *)malloc(ng*sizeof(float));
		
		for(i=0; i< ng; i++) {
			e_ave[i]= e_ave_show[i]= 0.0;
		}

		f_array= (float   *)malloc(nfft*sizeof(float));
		cur_hist= (float  *)malloc(nfft*sizeof(float));
		phi_hist= (float **)malloc(2*sizeof(float *));
		phi_hist[0]= (float *)malloc(nfft*sizeof(float));
		phi_hist[1]= (float *)malloc(nfft*sizeof(float));
		Local_t_array = (float *) malloc(nfft*sizeof(float));
		cur_fft= (float *)malloc(nfft*sizeof(float));
		phi_fft= (float *)malloc(nfft*sizeof(float));
		mphi_fft=(float *)malloc(nfft*sizeof(float));

		for(i=0; i< nfft; i++) {
			f_array[i]= i*df;
			cur_fft[i] = phi_fft[i] = mphi_fft[i] = 0.0;
		}
	}   /* End of if (nfft) */

	/** Allocate particle arrays **/
	x  = (double **)malloc(nsp*sizeof(double *));      /* double precision! */
	vx = (float **)malloc(nsp*sizeof(float *));
	vy = (float **)malloc(nsp*sizeof(float *));
	vz = (float **)malloc(nsp*sizeof(float *));
	c_omt = (int **)malloc(nsp*sizeof(int *));         /* for i_t=2 */
	gridref_r = (float *)malloc(nsp*sizeof(float));	/* reflection on RIGHT grid border */
	gridref_l = (float *)malloc(nsp*sizeof(float));	/* reflection on LEFT grid border */

	if (slow || dwindow) {
		xd = (float **)malloc(nsp*sizeof(float *));        /* For Xwindows */
	}

	for(i=0; i<nsp; i++) {
		x[i]  = (double *)malloc(MAXLEN*sizeof(double));
		vx[i] = (float *)malloc(MAXLEN*sizeof(float));
		vy[i] = (float *)malloc(MAXLEN*sizeof(float));
		vz[i] = (float *)malloc(MAXLEN*sizeof(float));
		c_omt[i] = (int *)malloc(MAXLEN*sizeof(int));
		
		if (slow || dwindow) {
			xd[i] = (float *)malloc(MAXLEN*sizeof(float));     /* For Xwindows */
		}

		for (k=0; k<MAXLEN; k++) {
			c_omt[i][k] = 0;
		}
	}
	
	if (nsp && (!x[nsp-1] || !vx[nsp-1] || !vy[nsp-1] || !vz[nsp-1])) {
		puts("START: x, v malloc failed");
		exit(1);
	}
} /* End of ALSP */



/*******  Free space, allocated in the beginning of program  ********/

void freespace() {
	int i, isp, k, l;

	
	for(i=0; i<nsp; i++) {
		free(x[i]);
		free(vx[i]);
		free(vy[i]);
		free(vz[i]);
		free(c_omt[i]);
		
		if (slow || dwindow) {
			free(xd[i]);     /* For Xwindows */
		}
	}

	if (slow || dwindow) {
		free(xd);        /* For Xwindows */
	}
	
	/** free particle arrays **/
	free(x);      /* double precision! */
	free(vx);
	free(vy);
	free(vz);
	free(c_omt);         /* for i_t=2 */
	free(gridref_r);	/* reflection on RIGHT grid border */
	free(gridref_l);	/* reflection on LEFT grid border */

	/** Allocate space for frequency history diagnostics **/
	if(nfft) {
		free(f_array);
		free(cur_hist);
		free(phi_hist[0]);
		free(phi_hist[1]);
		free(phi_hist);
		free(Local_t_array);
		free(cur_fft);
		free(phi_fft);
		free(mphi_fft);

		free(e_ave);
		free(e_ave_show);
		free(phi_ave);
		free(phi_ave_show);

		for(isp=0; isp<nsp; isp++) {
			free(srho_ave[isp]);
			free(srho_ave_show[isp]);
			free(Vx_ave[isp]);
			free(Vx_ave_show[isp]);
			free(Vpar_ave[isp]);
			free(Vpar_ave_show[isp]);
			free(Tx_ave[isp]);
			free(Qh_ave[isp]);
			free(Qtot_ave[isp]);
			free(Tx_ave_show[isp]);
			free(Qh_ave_show[isp]);
			free(Qtot_ave_show[isp]);
			free(T_dif_ave[isp]);
			free(T_dif_ave_show[isp]);

		}
	
		free(srho_ave);
		free(srho_ave_show);
		free(Vx_ave);
		free(Vx_ave_show);
		free(Vpar_ave);
		free(Vpar_ave_show);
		free(Tx_ave);
		free(Qh_ave);
		free(Qtot_ave);
		free(Tx_ave_show);
		free(Qh_ave_show);
		free(Qtot_ave_show);
		free(T_dif_ave);
		free(T_dif_ave_show);

	}   /* End of if (nfft) */

	if (slow || stopfl) {
		for(i=0; i<nsp; i++) {
			free(kes_hist[i]);
			free(jwall_hist[i]);
			free(Tx_hist[i]);
			free(T_dif_hist[i]);
			free(Vxt[i]);
		}

		free(ese_hist);
		free(Teq);
		free(kes_hist);
		free(jwall_hist);
		free(Tx_hist);
		free(T_dif_hist);
		free(Vxt);
	}
	
	for(i=0; i<nsp; i++) {
		if (slow || stopfl || Ef_flag) {
			free(Eflhs_hist[i]);
			free(Efrhs_hist[i]);
			free(Tlhs_hist[i]);
			free(Trhs_hist[i]);
		}

		free(np_hist[i]);
		free(F_hist[i]);
		free(F_histRHS[i]); 
		free(gridref_r_hist[i]);
		free(gridref_l_hist[i]);
	}
	
	if (slow || stopfl || Ef_flag) {
		for (i=0; i < Nnum; i++) {
			free(Halhs_hist[i]);
			free(Harhs_hist[i]);
		}

		free(Eflhs_hist);
		free(Efrhs_hist);
		free(Tlhs_hist);
		free(Trhs_hist);
		free(Halhs);
		free(Harhs);
		free(Halhs_temp);
		free(Harhs_temp);
		free(Halhs_hist);
		free(Harhs_hist);
	}

	free(t_array);
	free(np_hist);
	free(com_cur_hist);
	free(com_phi_hist[0]);
	free(com_phi_hist[1]);
	free(com_phi_hist);
	free(F_hist);
	free(F_histRHS);
	free(gridref_r_hist);
	free(gridref_l_hist);

	if (slow || stopfl) {
		/** Free space for diagnostics **/
		if (icollns) {
			for (isp=0; isp < Nnum; isp++) {
				free(chrgxrate[isp]);
			}

			free(chrgxrate);
		}

		/** Free space for collision diagnostics **/
		if (ecollisional) {
			for (isp=0; isp< Nnum; isp++) {
				free(ionrate[isp]);
				free(extrate[isp]);
			}

			free(ionrate);
			free(extrate);
		}

		for(isp=0; isp<nsp; isp++) {
			free(Vx[isp]);
			free(Vpar[isp]);
			free(Vy[isp]);
			free(Vz[isp]);
			free(Tx[isp]);
			free(Qh[isp]);
			free(Qtot[isp]);
			free(T_dif[isp]);
		}

		free(emin);
		free(de);
		free(dtheta);
		free(Eflhs);
		free(Eflhs_temp);
		free(Efrhs);
		free(Efrhs_temp);
		free(Tlhs);
		free(Tlhs_temp);
		free(Trhs);
		free(Trhs_temp);
		free(de_fvx);
		free(vx_min);

		free(nbin);
		free(nbin_fvx);
		
		free(Vx);
		free(Vpar);
		free(Vy);
		free(Vz);
		free(Tx);
		free(Qh);
		free(Qtot);
		free(T_dif);
		free(Qscale);
		free(Qwscale);

	/* end of (slow || stopfl) */
	} else {
		if (qtp>2) {          /* for effective ionization and/or movie */
			for(isp=0; isp<nsp; isp++) {
				free(Tx[isp]);
				
				if (qtp==3)	{
					free(Vpar[isp]);
					free(Qtot[isp]);
					free(Qh[isp]);
				}
			}

			if (qtp==3) {
				free(Vpar);
				free(Qtot);
				free(Qh);
				free(Qscale);
			}

			free(Vxef);
			free(Vyef);
			free(Vzef);
			free(Tx);
		}

		if (Ef_flag) {        /* for energy flux */
			free(Eflhs);
			free(Eflhs_temp);
			free(Efrhs);
			free(Efrhs_temp);
			free(Tlhs);
			free(Tlhs_temp);
			free(Trhs);
			free(Trhs_temp);
			free(Qwscale);
		}
	}   /* end of else */

	/* Surface physics */
	for(isp=0; isp<nsp; isp++) {
		for(i=0; i<MAXSRF; i++) {
			free(injtype[isp][i]);
			free(seec[isp][i]);
			free(Tseec[isp][i]);
			free(Ew[isp][i]);
			free(E0[isp][i]);
		}
	}

	for(isp=0; isp<nsp; isp++) {
		free(srf[isp]);
		free(secsp[isp]);
		free(sectype[isp]);

		free(injtype[isp]);
		free(seec[isp]);
		free(Tseec[isp]);
		free(Ew[isp]);
		free(E0[isp]);
	}

	free(srf);
	free(secsp);
	free(sectype);
	free(injtype);
	free(seec);
	free(Tseec);
	free(Ew);
	free(E0);

	
	
	if (Nnum) {
		/* Charged-neutral particle collisions */
		k= 1;
		
		for (isp=0; isp<Nnum; isp++) {
			for (i=0; i < Nnum; i++) {
				if (icollisional[isp][i] && k) {
					for (l=0; l < Nnum; l++) {
						free(nimiu[l]);
						free(sctang[l]);
					}

					free(nimiu);
					free(sctang);
					k= 0;
				}
			}
		}

		if ((ecollisional) && (Nnum)) {
			for (i=0; i<3; i++) {
				free(sctane[i]);
			}

			free(sctane); // For scat. angle
			free(Nexel); // Max. number of excit. coll.
		}

		
		for (i=0; i < Nnum; i++) {
			free(Nneut[i]);
			free(Tneut[i]);

			if (slow) {
				free(Tneutshow[i]);
			}
		}

		free(Nneut);
		free(Tneut);

		if (slow) {
			free(Tneutshow);
		}

	}  /* end of charged-neutral particle collisions */
	

	if (ccc) {              /* Coulomb collisions */
		for(isp=0; isp<=ccc; isp++) {
			for (i=0; i< nc; i++) {
				free(Ldav[isp][i]);
			}

			free(Ldav[isp]);
		}

		free(Ldav);
		free(pp);
		free(LL);
		free(Z);
		free(miu);
	}


	if (bmag > 0.0)	{
		free(tkv);
		free(tx);
		free(tz);
		free(sx);
		free(sz);
		free(count_max);
	}

	
	if (Ey_flag) {               /* Field parallel to the wall */
		free(Ey);
		free(E1y);
	}

	if (field) {
		for (isp=0; isp < nsp; isp++) {
			free(av_e[isp]);
		}

		free(av_e);
	}

	if (sourcD1) {           /* Second particle source */
		if (ht1>0.) {              /* "Heating" of the source */
			free(heat1);
		}

		free(vspart1);
		free(vspar01);
		free(vsperpt1);
		free(vsperp01);
		free(Src_con1);
		free(Stemp1);
	}
	
	if (sourcD) {                /* Particle source */
		free(vspart);
		free(vspar0);
		free(vsperpt);
		free(vsperp0);
		free(Src_con);
		free(Stemp);
		
		if (ht>0.) {               /* "Heating" of the source */
			free(heat);
		}
		
		if( sourcD == 4 ){
			/* allocate additional arrays... */
			free(totalLen) ;
			free(kLen) ;
			free(phi_six);
		}
		
	} else if (ccc) {  /* For calculation of Landau log. */
		free(vspart);
	}

	/* Initial loading */
	free(vpart);
	free(vpar0);
	free(vperpt);
	free(vperp0);

	for (isp=0; isp < nsp; isp++) 	{
		free(enter[isp]);
		free(vnpart[isp]);
		free(vnperpt[isp]);
		free(vnpar0[isp]);
		free(vnperp0[isp]);
		free(parvell[isp]);
		free(parvelr[isp]);
	}

	/* Injection from walls */
	free(enter);
	free(vnpart);
	free(vnperpt);
	free(vnpar0);
	free(vnperp0);
	free(parvell);
	free(parvelr);


	for(isp=0; isp<nsp; isp++) {
		free(srho[isp]);
	}

	free(Nsec);
	free(srho);

	free(gridref_r);		/* current reflection on RIGHT grid border */
	free(gridref_r0);		/* initial reflection on RIGHT grid border */
	free(gridref_rsramp);	/* size step of reflection on RIGHT grid border - values is multiplied by time */
	free(gridref_l);		/* current reflection on LEFT grid border */
	free(gridref_l0);		/* initial reflection on LEFT grid border */
	free(gridref_lsramp);	/* size step of reflection on LEFT grid border - values is multiplied by time */
	
	/** free space for general parameters **/
	free(phi);
	free(rho);
	free(e);

	free(fr_E);
	free(jwall);
	free(qm);
	free(q);
	free(m);
	free(Escale);
	free(Tscale);
	free(Fscale);
	free(np);
	free(nstep);
	free(wt);
	free(ht_sp);
	free(dinj);
	free(i_t);
	free(inject);
	free(F_temp);
	free(F_tempRHS);
	free(FLHS);
	free(FRHS);

	free(x_array);

}


/***************************************************************/
