import BoxTools.Util;
public class Operator{
	static final int SPACE_DIM=Util.SPACE_DIM;
	final static int highBoundary=Util.HIGH_SIDE;
	final static int lowBoundary=Util.LOW_SIDE;
	RectDomain<SPACE_DIM> m_domain;
	RectDomain<SPACE_DIM> [] local m_boundaries=new RectDomain<SPACE_DIM> [2*SPACE_DIM];
	public Operator(RectDomain<SPACE_DIM> a_domain){
		m_domain=a_domain;
		int dir, side,i;
		for (i=0;i<2*SPACE_DIM;i++){
			dir=(i+2)/2;
			side=(i+2)%2;
			if (side==highBoundary)
				m_boundaries[i]=m_domain.border(1,dir,1);
			else
				m_boundaries[i]=m_domain.border(1,-dir,1);	
		}

	}
	public final local void applyOp(template BoxedArray<double> local a_Arho, template BoxedArray<double> local a_rho){
		int dir, side,i;
		Point<SPACE_DIM> [1d] local disp =new Point<SPACE_DIM> [[1:SPACE_DIM]];
		for (i=1;i<=SPACE_DIM;i++) disp[i]=Point<SPACE_DIM>.direction(i);
		double [SPACE_DIM d] local TiArray1, TiArray2;
		TiArray1=a_rho.getLocalArray();
		TiArray2=a_Arho.getLocalArray();

		for (i=0;i<2*SPACE_DIM;i++){
			dir=(i+2)/2;
			side=(i+2)%2;
			if (side==highBoundary)
				foreach (point in  a_boundaries[i]) TiArray1[point]=-TiArray1[point-disp[dir]];
			else
				foreach (point in  a_boundaries[i]) TiArray1[point]=-TiArray1[point+disp[dir]];	
		}
		RectDomain<SPACE_DIM> box=a_Arho.domain();
		switch (SPACE_DIM){
				   case 1:
				     foreach (point in box)
				       TiArray2[point]=(TiArray1[point+disp[1]]+TiArray1[point-disp[1]]-2*TiArray1[point])/h2;
				     break;
				   case 2:
				     foreach (point in box)
				       TiArray2[point]=(TiArray1[point+disp[1]]+TiArray1[point-disp[1]]+TiArray1[point+disp[2]]+
						     TiArray1[point-disp[2]]-4*TiArray1[point])/h2;
				     break;
				   case 3:
				     foreach (point in box)
				       TiArray2[point]=(TiArray1[point+disp[1]]+TiArray1[point-disp[1]]+TiArray1[point+disp[2]]+
						     TiArray1[point-disp[2]]+TiArray1[point+disp[3]]+TiArray1[point-disp[3]]-
						     6*TiArray1[point])/h2;
				 }
	return;
	}
}