Hi Vijay,

Thanks for your first comments.

Here is a first version of arrays on 3D distributions including ghost regions.

I started from the MortonDist class of Josh Milthorpe et al.
Basically, I defined a distribution by 3D blocks (if I'm not mistaken, X10 only provides 1D and 2D block distributions).
Each local block is extended by a few planes on each face.

The distribution defines the specific methods
- getLocal (block without ghost regions) and
- getExtended (block with supplementary planes as ghost regions)

I only tried it on a simple example (see below).

I would be very interested if you or any other x10 expert can spent some time to send me comments or suggestions.

Thanks in advance,
Marc

Yes this should be straightforward to do -- provided you are willing to
write a wrapper clause.

In fact this would be a good exercise to work through.

You may want to think carefully about the methods on the class. e.g. it
might be best to view the new class as having some synchronous methods
(called simultaneously on multiple threads, one at each place) for the
communication of ghost regions. Once this is done, the remote state will
be local and can now be operated upon by purely local code.

Asynchronous boundary exchangers are also possible.

However, not sure why you would like the new class to subclass Array or
Distr
On 1/20/2012 6:35 AM, TAJCHMAN Marc wrote:
>  Hi,
>
>  I have the following problem
>
>  Let's consider a distributed array : a(0) .. a(n)  (1D version for 
simplicity):
>
>  a(1) .. a(10) on place 0
>  a(11) .. a(20) on place 1
>  a(21) .. a(30) on place 2
>  ..
>
>  I want to apply a stencil (e.g. from the heat equation : b_(i) = f(a_(i+1), 
a(i), a_(i-1))).
>  where b has the same distribution as a
>
>  Is it possible to define a class (that inherits from array) so that
>  __________________________________________________________
>  - at each place p, the local part of the array a(1 + 10*p) ... a(10*(p+1)) 
is surrounded
>     by a layer of "ghost components" that are copies of a(10*p) and 
a(10*(p+1) + 1)
>
>  __________________________________________________________
>  - at a place p:
>
>  if you write :
>
>           a(1+10*p) = a(10*p),
>
>  you read the local copy (ghost value) of the distant component a(10*p)
>
>    if you write :
>
>            a(1+10*p) = at (p-1) a(10*p),
>
>  you have a remote access to the real component
>
>  __________________________________________________________
>  - the class has a method (e.g. ghosts_synchronize) that efficiently
>     duplicates selected distant values in ghost values
>
>  Maybe this kind of arrays have already been implemented somewhere ?
>
>  __________________________________________________________
>  Thanks in advance,
>  Marc


--
__________________________________________________________
Marc Tajchman
Laboratoire de Génie Logiciel et de Simulation
Software Engineering and Computer Simulation Laboratory
CEA-DEN/DANS/DM2S/STMF/LGLS
91191 Gif-sur-Yvette cedex

Tél : +33/0 1 69 08 73 27
Fax : +33/0 1 69 08 96 96
E-mail:/marc.tajch...@cea.fr <mailto:marc.tajch...@cea.fr>/
__________________________________________________________

/Ce message électronique et tous les fichiers attachés qu'il contient sont confidentiels et destinés exclusivement à l'usage de la personne à laquelle ils sont adressés. Si vous avez reçu ce message par erreur, merci d'en avertir immédiatement son émetteur et de ne pas en conserver de copie./

/This e-mail and any files transmitted with it are confidential and intended solely for the use of the individual to whom they are addressed. If you have received this e-mail in error, please inform the sender immediately, without keeping any copy thereof./
/*
  3D Block distribution with ghost regions on the border of each block.
  A ghost region is a local duplicate of points located on surrounding 
  (distant) blocks, neighbors of local points.
*/
public class Block3DDist extends Dist(3) {

  /* Group of places where the distribrution lives */
  private val pg:PlaceGroup;

  /* Number of subdomains in each direction */
  private val nDom: Array[int](1);

  /* Size of a extended local region (including ghost regions) */
  private transient var extendedRegion:Region(3);

  /* Standard local sizes (in each direction) of a subdomain (except for the
     last subdomain in each direction) */
  private val stdSize: Array[int](1);

  /* Depth of ghost regions in each direction */
  private val overlay : int;

  /* Factories */
  public static def makeBlock3D(R:Region(3), nDom: Array[int](1), overlay:int)
  {
    return new Block3DDist(R, nDom, overlay, PlaceGroup.WORLD);
  }

  public static def makeBlock3D(R:Region(3), nDom:Array[int](1), overlay:int, 
                                pg:PlaceGroup)
  {
    return new Block3DDist(R, nDom, overlay, pg);
  }

  /* Class ctor : computes the standard size of a subdomain (excluding
     ghost region */
  def this(R:Region(3), nD:Array[int](1), overlay:int, pg:PlaceGroup) 
    : Block3DDist
  {
    super(R);

    if (pg.numPlaces() != (nD(0) * nD(1) * nD(2)))
      throw new RuntimeException("wrong number of subdomains" + nD);

    this.pg = pg;
    this.nDom = new Array[int](rank, (i:int) => nD(i));
    var _Size : Array[int](1) = new Array[int]
      (rank, (i : Int) => 
       Math.floor((R.max(i)-R.min(i)+1)/nD(i)) as Int);
    for (i in 0..(rank-1)) {
      if (_Size(i) * nD(i) < R.max(i) - R.min(i) + 1)
        _Size(i) += 1;
    }
    this.stdSize = _Size;
    this.overlay = overlay;
  }

  public def places():PlaceGroup = pg;
  
  public def numPlaces():Int = pg.numPlaces();
 
  /* Given the (id,jd,kd) index of a subdomain, return the place 
     where it lives */
  private def compose(val idx:Array[int](1)) : Place {

    if ((idx(0) < 0) || (idx(0) >= nDom(0)) ||
        (idx(1) < 0) || (idx(1) >= nDom(1)) ||
        (idx(2) < 0) || (idx(2) >= nDom(2))) 
      throw new RuntimeException("wrong subdomain index" + idx);

    val q = idx(0) + nDom(0)*(idx(1) + nDom(1)*idx(2));
    // Console.ERR.println("compose " + idx + " = " + q);
    return pg(q);
  }

  /* From a place, computes the (id,jd,kd) index of the subdomain 
     located at this place 

     compose and decompose are inverse operations:
          decompose(compose((id,jd,kd)) = (id,jd,kd)
          compose(decompose(place)) = place
  */
  private def decompose(val p:Place) {
    
    val id = pg.indexOf(p);
    val idx = new Array[int]((0..2));
    
    var i:int = id;
    
    idx(2) = i/(nDom(0)*nDom(1));
    i -= idx(2) * nDom(0)*nDom(1);
    
    idx(1) = i/nDom(0);
    i -= idx(1) * nDom(0);
    
    idx(0) = i;
    return idx;
  }

  /* Returns the extended local region for a place */
  private def regionForPlace(place:Place) : Region{self.rank==this.rank} {

    val idxDom = decompose(place);
      
    val newMin = new Array[Int](rank, (i:Int) => 
                                region.min(i) + stdSize(i)*idxDom(i)-overlay);
    val newMax = new Array[Int](rank, (i:Int) =>
                                idxDom(i) < nDom(i) - 1 
                                ? region.min(i) + stdSize(i)*(idxDom(i)+overlay)
                                : region.max(i) + overlay);
   
    return Region.makeRectangular(newMin, newMax);
  }

  public def regions():Sequence[Region(rank)] {
    return new Array[Region(rank)]
      (pg.numPlaces(), (i:int)=>regionForPlace(pg(i))).sequence();
  }

  /* Returns the extended local region for a place */
  public def get(p:Place): Region(rank) {
    if (p == here) {
      if (extendedRegion == null) {
        extendedRegion = regionForPlace(here);
      }
      return extendedRegion;
    } else {
      return regionForPlace(p);
    }
  }

  /* Returns the local region for a place (excluding ghost regions) */
  public def getLocal(p:Place):Region(rank) {
    val R = get(p);
    return (R.min(0) + overlay) .. (R.max(0) - overlay) *
      (R.min(1) + overlay) .. (R.max(1) - overlay) *
      (R.min(2) + overlay) .. (R.max(2) - overlay);
  }

  /* Returns the extended region for a place (including ghost regions) */
  public def getExtended(p:Place):Region(rank) {
    return get(p);
  }

  /* Returns the place where a point lives in the distribution */
  public operator this(pt:Point(rank)):Place {

    /* _id_p : (id,jd,kd) index of the subdomain containing the point */
    val _id_p = new Array[Int]
      (rank, (i : Int) => 
       Math.floor((pt(i)-region.min(i))/stdSize(i)) as Int);
    return compose(_id_p);
  }


  public def restriction(r:Region(rank)):Dist(rank) {
    throw new UnsupportedOperationException("restriction(r:Region(rank))");
  }

  public def restriction(p:Place):Dist(rank) {
    return Dist.makeConstant(this.get(p));
  }

  public def offset(pt:Point(rank)):Int {
    return get(here).indexOf(pt);
  }

  public def maxOffset() {
    val r = get(here);
    return r.size()-1;
  }

  /* Places containing neighbor subdomains of the subdomain in place p.

     If the local subdomain is on a border of the global region in some
     directions, neighbor subdomain is the local subdomain itself for these
     directions

     Q: Cannot return a null place ? */
  public def neighbors(p:Place): Sequence[Place]
  {
    val idx = decompose(p);
    val idxx = new Array[Int](rank, (i : Int) => idx(i));
    val n = new Array[Place](0..(2*rank-1));

    for ([i] in 0..(rank-1)) {
      idxx(i) = idx(i) - 1;
      n(2*i)   = (idxx(i) >= 0)      ? compose(idxx) : p;
      idxx(i) = idx(i) + 1;
      n(2*i+1) = (idxx(i) < nDom(i)) ? compose(idxx) : p;
      idxx(i) = idx(i);
    }
    return n.sequence();
  }
 
  /* Printable representation of the distribution */
  public def toString(): String {
    var s: String = "Block3DDist(\n";
    var first: boolean = true;
    for (p:Place in places()) {
      val r = getLocal(p);
      val r_ext = getExtended(p);
      if (!r.isEmpty()) {
        s += "    ";
        if (!first) s += ",";
        s += "" + r + "(ext. " + r_ext + ") ->" + p.id + "\n";
        first = false;
      }
    }
    s += ")";
    return s;
  }

  /* For a distributed array defined on this distribution, 
     generic function that copy border values in each subdomain
     tho ghost regions of surrounding subdomains */
  static public def synchronize[T](A : DistArray[T](3))
  {
    val d = A.dist as Block3DDist;
        
    finish
      for (p in d.places()) at (p) async {
          val n = d.neighbors(p);
          val r = d.getLocal(p);
          val r_exchange = new Array[Region(d.rank)](0..(2*d.rank));
          r_exchange(0) 
            = (r.min(0) .. r.min(0))
            * (r.min(1) .. r.max(1))
            * (r.min(2) .. r.max(2));
          r_exchange(1) 
            = (r.max(0) .. r.max(0))
            * (r.min(1) .. r.max(1))
            * (r.min(2) .. r.max(2));
          r_exchange(2) 
            = (r.min(0) .. r.max(0))
            * (r.min(1) .. r.min(1))
            * (r.min(2) .. r.max(2));
          r_exchange(3) 
            = (r.min(0) .. r.max(0))
            * (r.max(1) .. r.max(1))
            * (r.min(2) .. r.max(2));
          r_exchange(4) 
            = (r.min(0) .. r.max(0))
            * (r.min(1) .. r.max(1))
            * (r.min(2) .. r.min(2));
          r_exchange(5) 
            = (r.min(0) .. r.max(0))
            * (r.min(1) .. r.max(1))
            * (r.max(2) .. r.max(2));
                    
          for ([l] in 0..(2*d.rank-1)) async {
              if (n(l) != p) {
                val q = n(l);
                val r_ex = r_exchange(l);
                val A_local = new Array[T]
                  (r_ex, ([i,j,k]:Point) => A(i,j,k));
                val A_remote = at (q) new RemoteArray(A_local);
                finish Array.asyncCopy(A_local, A_remote);
                at (q) {
                  for ([i,j,k] in r_ex)
                    A(i,j,k) = A_local(i,j,k);
                }
              }
            }
                    
        }
  }

}
import Block3DDist;
import x10.io.Printer;
import x10.io.File;

public class main_TestDist
{

  public static def output(message : String, W: DistArray[double](3))
  {
    val r = W.region;
    val dx = 1.0/(r.max(0) - r.min(0) + 1);
    val dy = 1.0/(r.max(1) - r.min(1) + 1);
    val dz = 1.0/(r.max(2) - r.min(2) + 1);

    val fileName = "out_" + Place.MAX_PLACES.toString() 
      + "_" + message + ".vts";
    val out = new File(fileName).printer();
    out.println("<?xml version=\"1.0\"?>");
    out.println("<VTKFile type=\"StructuredGrid\">");
    out.println("<StructuredGrid WholeExtent=\"" 
                + " " + (r.min(0)-1) + " " + r.max(0)  
                + " " + (r.min(1)-1) + " " + r.max(1) 
                + " " + (r.min(2)-1) + " " + r.max(2)
                + "\">");
    out.println("<Piece Extent=\""
                + " " + (r.min(0)-1) + " " + r.max(0) 
                + " " + (r.min(1)-1) + " " + r.max(1) 
                + " " + (r.min(2)-1) + " " + r.max(2)
                + "\">");

    out.println("<Points>");
    out.println("<DataArray type=\"Float32\" " + 
                "format=\"ascii\" NumberOfComponents=\"3\">");

    for ([k] in r.min(2)..(r.max(2)+1))
      for ([j] in r.min(1)..(r.max(1)+1))
        for ([i] in r.min(0)..(r.max(0)+1)) {
          val x = (i-r.min(0))*dx, 
            y = (j-r.min(1))*dy, 
            z = (k-r.min(2))*dz;
          out.println(" " + x + " " + y + " " + z);
        }
    
    out.println("</DataArray>");
    out.println("</Points>");

    out.println("<CellData Scalars=\"" + message +"\">\n");
    out.println("<DataArray type=\"Float32\" Name=\""
                + message + "\" format=\"ascii\">");
    for ([k] in r.min(2)..r.max(2))
      for ([j] in r.min(1)..r.max(1))
        for ([i] in r.min(0)..r.max(0)) {
          val v = at (W.dist(i,j,k)) W(i,j,k);
          out.print(" " + v);
        }
    out.println("</DataArray>");
    
    out.println("</CellData>");
    out.println("</Piece>");
    out.println("</StructuredGrid>");
    out.println("</VTKFile>");
    out.close();
  }

  public static def iterate(U: DistArray[double](3), V: DistArray[double](3))
  {
    finish
      for (p in U.dist.places()) 
        async at (p) {
          val r = (U.dist as Block3DDist).getLocal(p);
          for ([i,j,k] in r)
            V(i,j,k) = (U(i+1,j,k) + U(i-1,j,k) + 
                        U(i,j+1,k) + U(i,j-1,k) +
                        U(i,j,k+1) + U(i,j,k-1))/6.0;
        }
    Block3DDist.synchronize(V);
  }
 
  public static def main(args:Array[String](1)) {

    val start_time = System.nanoTime();
    val n : Int;
    if (args.size > 0)
      n = Int.parse(args(0));
    else
      n = 10;

    val nDom = new Array[int](0..2);
    for (i in 0..2)
      if (args.size > i+1)
        nDom(i) = Int.parse(args(i+1));
      else
        nDom(i) = 1;

    Console.ERR.println("nDomains = " + 
                        nDom(0) + " x " + nDom(1) + " x " + nDom(2));

    val R = 1..n * 1..n * 1..n;
    Console.ERR.println("Global region = " + R);

    val D = Block3DDist.makeBlock3D(R, nDom, 1);
    Console.ERR.println("Distribution D = " + D + " D.region" + D.region);

    for (p in Place.places())
      Console.ERR.println("D(" + p + ") = " + D(p) + " local " + D.getLocal(p));

    val start_time_init = System.nanoTime();
    val U = DistArray.make[double](D);
    val V = DistArray.make[double](D);

    finish
      for (p in Place.places()) 
        at (p) async {
          val r = (U.dist as Block3DDist).getLocal(p);
          for ([i,j,k] in r) {
            U(i,j,k) = ((Math.abs(i-n/2) < n/4) && 
                        (Math.abs(j-n/2) < n/4) && 
                        (Math.abs(k-n/2) < n/4)) 
              ? 1.0 : 0.0;
          }
        }
    Block3DDist.synchronize(U);
    Console.ERR.println("Time(init)      = " 
                        + (1e-9*(System.nanoTime() - start_time_init)));

    // output("U0", U);   

    val start_time_iterate = System.nanoTime();
    var b:Boolean = true;
    
    Console.ERR.print("iteration ");
    for ([it] in 0..100) {
      Console.ERR.print(" " + it);
      if (b) {
        iterate(U,V);
        //            output("V", V);   
      }
      else {
        iterate(V,U);
        //            output("U", U);   
      }
      b = !b;
    }
    Console.ERR.println();
    Console.ERR.println("Time(iterate)   = " 
                        + (1e-9*(System.nanoTime() - start_time_iterate)));

    Console.ERR.println("Time(total)     = " 
                        + (1e-9*(System.nanoTime() - start_time_init)));

    output("U", b ? U : V);
  }
}
------------------------------------------------------------------------------
Keep Your Developer Skills Current with LearnDevNow!
The most comprehensive online learning library for Microsoft developers
is just $99.99! Visual Studio, SharePoint, SQL - plus HTML5, CSS3, MVC3,
Metro Style Apps, more. Free future releases when you subscribe now!
http://p.sf.net/sfu/learndevnow-d2d
_______________________________________________
X10-users mailing list
X10-users@lists.sourceforge.net
https://lists.sourceforge.net/lists/listinfo/x10-users

Reply via email to