freepooma-devel
[Top][All Lists]
Advanced

[Date Prev][Date Next][Thread Prev][Thread Next][Date Index][Thread Index]

[PATCH] Add PETSc wrapper.


From: Richard Guenther
Subject: [PATCH] Add PETSc wrapper.
Date: Tue, 23 Mar 2004 10:30:20 +0100 (CET)

Hi!

This adds support for using the PETSc library of sparse (non-)linear
(MPI-)parallel solvers with POOMA.  It does so by excercising the PETSc
distributed array (DA) API for interfacing with POOMA engines.

Accompanied with an example in Solvers/PETSc.  Tested with serial
PETSc/POOMA and MPI PETSc/POOMA.

Ok?

Richard.


2004Mar23  Richard Guenther <address@hidden>

        * src/Transform/PETSc.h: new.
        examples/Solvers/PETSc/PETSc.cpp: new.
        examples/Solvers/PETSc/include.mk: new.
        examples/Solvers/PETSc/makefile: new.

diff -Nru a/r2/examples/Solvers/PETSc/PETSc.cpp 
b/r2/examples/Solvers/PETSc/PETSc.cpp
--- /dev/null   Wed Dec 31 16:00:00 1969
+++ b/r2/examples/Solvers/PETSc/PETSc.cpp       Tue Mar 23 10:25:34 2004
@@ -0,0 +1,100 @@
+// Example on how to use PETSc linear solvers from inside a Pooma program.
+// This solves the Poisson equation for a density distribution in a 3D
+// Array to create the corresponding gravitational potential.
+
+#include <iostream>
+#include "Transform/PETSc.h"
+
+#include "petscda.h"
+#include "petscksp.h"
+
+void doit()
+{
+  typedef Array<2, double, MultiPatch<GridTag, Remote<Brick> > > Array_t;
+  typedef Array_t::Layout_t::const_iterator layout_iter;
+
+  Interval<2> domain(8, 8);
+  Loc<2> blocks = makeRBlocks(domain, Pooma::contexts());
+  GridLayout<2> layout(domain, blocks, GuardLayers<2>(1), DistributedTag());
+
+  // Create DA
+  Pooma::PoomaDA da(layout, DA_NONPERIODIC, DA_STENCIL_STAR, 1);
+
+  // Assemble matrix
+  Mat jac;
+  DAGetMatrix(da, MATMPIAIJ, &jac);
+  {
+    int mx, my;
+    PetscScalar Hx, Hy, HxdHy, HydHx;
+    DAGetInfo(da,0,&mx,&my,PETSC_NULL,0,0,0,0,0,0,0);
+    Hx = 1.0 / (PetscReal)(mx-1); Hy = 1.0 / (PetscReal)(my-1);
+    HxdHy = Hx/Hy; HydHx = Hy/Hx;
+    int xs, ys, xm, ym;
+    MatStencil row,col[5];
+    DAGetCorners(da,&xs,&ys,PETSC_NULL,&xm,&ym,PETSC_NULL);
+    PetscScalar v[5];
+    for (int j=ys; j<ys+ym; j++){
+      for(int i=xs; i<xs+xm; i++){
+       row.i = i; row.j = j;
+       if (i==0 || j==0 || i==mx-1 || j==my-1) {
+         v[0] = 2.0*(HxdHy + HydHx);
+         MatSetValuesStencil(jac,1,&row,1,&row,v,INSERT_VALUES);
+       } else {
+         v[0] = -HxdHy;              col[0].i = i;   col[0].j = j-1;
+         v[1] = -HydHx;              col[1].i = i-1; col[1].j = j;
+         v[2] = 2.0*(HxdHy + HydHx); col[2].i = i;   col[2].j = j;
+         v[3] = -HydHx;              col[3].i = i+1; col[3].j = j;
+         v[4] = -HxdHy;              col[4].i = i;   col[4].j = j+1;
+         MatSetValuesStencil(jac,1,&row,5,col,v,INSERT_VALUES);
+       }
+      }
+    }
+  }
+  MatAssemblyBegin(jac,MAT_FINAL_ASSEMBLY);
+  MatAssemblyEnd(jac,MAT_FINAL_ASSEMBLY);
+
+  // problem
+  Array_t rh(layout), ph(layout);
+  rh(rh.totalDomain()) = 1.0;
+  ph(ph.totalDomain()) = 0.0;
+  Pooma::pinfo << "Density distribution:\n" << rh << std::endl;
+
+  // assemble rhs and solution
+  Vec pph, prh;
+  DACreateGlobalVector(da, &prh);
+  VecDuplicate(prh, &pph);
+
+  da.assign(prh, rh);
+
+  // create solver
+  KSP ksp;
+  KSPCreate(PETSC_COMM_WORLD, &ksp);
+  KSPSetTolerances(ksp, 1e-10, 1e-5, 1e5, 10);
+  KSPSetOperators(ksp, jac, jac, DIFFERENT_NONZERO_PATTERN);
+  KSPSetFromOptions(ksp);
+
+  // solve and copy solution
+  KSPSetRhs(ksp, prh);
+  KSPSetSolution(ksp, pph);
+  KSPSolve(ksp);
+  da.assign(ph, pph);
+
+  // show soultion
+  Pooma::pinfo << "Solution:\n" << ph << std::endl;
+}
+
+int main(int argc, char **argv)
+{
+  Pooma::initialize(argc, argv);
+
+  PetscSetCommWorld(MPI_COMM_WORLD);
+  PetscInitialize(&argc, &argv, NULL, NULL);
+
+  doit();
+
+  // cleanup
+  PetscFinalize();
+  Pooma::finalize();
+  return 0;
+}
+
diff -Nru a/r2/examples/Solvers/PETSc/include.mk 
b/r2/examples/Solvers/PETSc/include.mk
--- /dev/null   Wed Dec 31 16:00:00 1969
+++ b/r2/examples/Solvers/PETSc/include.mk      Tue Mar 23 10:25:34 2004
@@ -0,0 +1,59 @@
+# Generated by mm.pl: Mon Mar  9 13:58:39 MST 1998
+# ACL:license
+#  ----------------------------------------------------------------------
+#  This software and ancillary information (herein called "SOFTWARE")
+#  called POOMA (Parallel Object-Oriented Methods and Applications) is
+#  made available under the terms described here.  The SOFTWARE has been
+#  approved for release with associated LA-CC Number LA-CC-98-65.
+#
+#  Unless otherwise indicated, this SOFTWARE has been authored by an
+#  employee or employees of the University of California, operator of the
+#  Los Alamos National Laboratory under Contract No. W-7405-ENG-36 with
+#  the U.S. Department of Energy.  The U.S. Government has rights to use,
+#  reproduce, and distribute this SOFTWARE. The public may copy, distribute,
+#  prepare derivative works and publicly display this SOFTWARE without
+#  charge, provided that this Notice and any statement of authorship are
+#  reproduced on all copies.  Neither the Government nor the University
+#  makes any warranty, express or implied, or assumes any liability or
+#  responsibility for the use of this SOFTWARE.
+#
+#  If SOFTWARE is modified to produce derivative works, such modified
+#  SOFTWARE should be clearly marked, so as not to confuse it with the
+#  version available from LANL.
+#
+#  For more information about POOMA, send e-mail to address@hidden,
+#  or visit the POOMA web page at http://www.acl.lanl.gov/pooma/.
+#  ----------------------------------------------------------------------
+# ACL:license
+
+
+# Wrap make components from SHARED_ROOT and the current directory in the
+# proper order so that variables like ODIR have the correct directory-specific
+# value at the right moment.  See the included files for details of what they
+# are doing. This file should NOT be manually edited.
+
+# Set NEXTDIR, THISDIR and DIR_LIST
+include $(SHARED_ROOT)/include1.mk
+
+# Include list of subdirectories to process
+-include $(THISDIR)/subdir.mk
+
+# Set ODIR, PROJECT_INCLUDES, UNIQUE
+include $(SHARED_ROOT)/include2.mk
+
+# Set list of object files, relative to ODIR
+-include $(THISDIR)/objfile.mk
+
+# Set rules for the ODIR directory
+include $(SHARED_ROOT)/compilerules.mk
+
+# Remove current dir from DIR_LIST
+DIR_LIST :=$(filter-out $(firstword $(DIR_LIST)), $(DIR_LIST))
+
+
+# ACL:rcsinfo
+#  ----------------------------------------------------------------------
+#  $RCSfile: include.mk,v $   $Author: swhaney $
+#  $Revision: 1.3 $   $Date: 2000/03/07 13:15:37 $
+#  ----------------------------------------------------------------------
+# ACL:rcsinfo
diff -Nru a/r2/examples/Solvers/PETSc/makefile 
b/r2/examples/Solvers/PETSc/makefile
--- /dev/null   Wed Dec 31 16:00:00 1969
+++ b/r2/examples/Solvers/PETSc/makefile        Tue Mar 23 10:25:34 2004
@@ -0,0 +1,65 @@
+# Generated by mm.pl: Mon Mar  9 13:58:39 MST 1998
+# ACL:license
+#  ----------------------------------------------------------------------
+#  This software and ancillary information (herein called "SOFTWARE")
+#  called POOMA (Parallel Object-Oriented Methods and Applications) is
+#  made available under the terms described here.  The SOFTWARE has been
+#  approved for release with associated LA-CC Number LA-CC-98-65.
+#
+#  Unless otherwise indicated, this SOFTWARE has been authored by an
+#  employee or employees of the University of California, operator of the
+#  Los Alamos National Laboratory under Contract No. W-7405-ENG-36 with
+#  the U.S. Department of Energy.  The U.S. Government has rights to use,
+#  reproduce, and distribute this SOFTWARE. The public may copy, distribute,
+#  prepare derivative works and publicly display this SOFTWARE without
+#  charge, provided that this Notice and any statement of authorship are
+#  reproduced on all copies.  Neither the Government nor the University
+#  makes any warranty, express or implied, or assumes any liability or
+#  responsibility for the use of this SOFTWARE.
+#
+#  If SOFTWARE is modified to produce derivative works, such modified
+#  SOFTWARE should be clearly marked, so as not to confuse it with the
+#  version available from LANL.
+#
+#  For more information about POOMA, send e-mail to address@hidden,
+#  or visit the POOMA web page at http://www.acl.lanl.gov/pooma/.
+#  ----------------------------------------------------------------------
+# ACL:license
+
+# This file is user-editable
+
+PROJECT_ROOT = $(shell cd ../../..; pwd)
+include $(PROJECT_ROOT)/config/head.mk
+
+PASS=APP
+
+default:: PETSc
+
+.PHONY: PETSc
+
+PETSc:: $(ODIR)/PETSc
+
+include $(SHARED_ROOT)/tail.mk
+
+$(ODIR)/PETSc.o: PETSc.cpp
+       $(CXX) -o $@ -c $< $(RULE_CXXOPTS) $(RULE_INCLUDES) $(SUITE_DEFINES) \
+               -DPETSC_USE_EXTERN_CXX \
+               -I$(PETSC_DIR)/bmake/linux-gnu \
+               -I$(PETSC_DIR)/include
+
+$(ODIR)/PETSc: $(ODIR)/PETSc.o
+       $(CXX) -o $@ $^ $(RULE_LD_OPTS) \
+               -L$(PETSC_DIR)/lib/libg/linux-gnu \
+               -lpetscmat -lpetscksp -lpetscdm -lpetscvec -lpetscmat -lpetscdm 
-lpetsc -lpetscts -lpetscvec -lpetscsnes -lblas -llapack -lblas -lg2c -ldl 
-static
+
+# ACL:rcsinfo
+#  ----------------------------------------------------------------------
+#  $RCSfile: makefile,v $   $Author: julianc $
+#  $Revision: 1.6 $   $Date: 2000/07/21 20:43:14 $
+#  ----------------------------------------------------------------------
+# ACL:rcsinfo
diff -Nru a/r2/src/Transform/PETSc.h b/r2/src/Transform/PETSc.h
--- /dev/null   Wed Dec 31 16:00:00 1969
+++ b/r2/src/Transform/PETSc.h  Tue Mar 23 10:25:34 2004
@@ -0,0 +1,396 @@
+#ifndef POOMA_TRANSFORM_PETSC_H
+#define POOMA_TRANSFORM_PETSC_H
+
+// PETSc interfacing with POOMA using the PETSc DA interface.
+//
+// Copyright (c) 2004 by Richard Guenther <address@hidden>
+//
+// This file is in the public domain.
+
+/** @file
+ * @ingroup Utilities
+ * @brief
+ * Interfacing with the PETSc library of (non-)linear solvers.
+ *
+ * Interfacing supports the PETSc DA (distributed arrays) notion
+ * for creating (non-)linear solvers for implicit finite difference
+ * methods.  Using this wrappers you can fill your right-hand-side
+ * vector from a POOMA engine and transfer the result-vector to
+ * a POOMA engine.
+ *
+ * You are going to use the PetscDA class and its methods.
+ * See examples/Solver/PETSc for how to use this.
+ */
+
+#include "Pooma/Arrays.h"
+#include "petscda.h"
+
+
+template <class MeshTag, class T, class EngineTag>
+class Field;
+
+
+namespace Pooma {
+
+
+/**
+ * Helper to convert DALocalInfo domain info to appropriate
+ * Pooma Interval
+ */
+
+template <int Dim>
+struct PoomaDAGetDomain;
+
+template <>
+struct PoomaDAGetDomain<1> {
+  static inline
+  Interval<1> innerDomain(DALocalInfo& i)
+  {
+    return Interval<1>(i.xs, i.xs+i.xm-1);
+  }
+  static inline
+  Interval<1> totalDomain(DALocalInfo& i)
+  {
+    return Interval<1>(i.gxs, i.gxs+i.gxm-1);
+  }
+};
+
+template <>
+struct PoomaDAGetDomain<2> {
+  static inline
+  Interval<2> innerDomain(DALocalInfo& i)
+  {
+    return Interval<2>(Interval<1>(i.xs, i.xs+i.xm-1),
+                      Interval<1>(i.ys, i.ys+i.ym-1));
+  }
+  static inline
+  Interval<2> totalDomain(DALocalInfo& i)
+  {
+    return Interval<2>(Interval<1>(i.gxs, i.gxs+i.gxm-1),
+                      Interval<1>(i.gys, i.gys+i.gym-1));
+  }
+};
+
+template <>
+struct PoomaDAGetDomain<3> {
+  static inline
+  Interval<3> innerDomain(DALocalInfo& i)
+  {
+    return Interval<3>(Interval<1>(i.xs, i.xs+i.xm-1),
+                      Interval<1>(i.ys, i.ys+i.ym-1),
+                      Interval<1>(i.zs, i.zs+i.zm-1));
+  }
+  static inline
+  Interval<3> totalDomain(DALocalInfo& i)
+  {
+    return Interval<3>(Interval<1>(i.gxs, i.gxs+i.gxm-1),
+                      Interval<1>(i.gys, i.gys+i.gym-1),
+                      Interval<1>(i.gzs, i.gzs+i.gzm-1));
+  }
+};
+
+
+
+/**
+ * Helper to ease brick-engine -> vector copy
+ */
+
+template <int Dim>
+struct PoomaDACopy;
+
+template <>
+struct PoomaDACopy<1> {
+  template <class T>
+  static
+  void copy(Vec v, const Engine<1, T, Brick>& e)
+  {
+    PetscScalar *pa;
+    VecGetArray(v, &pa);
+    int idx=0;
+    Interval<1> d(e.domain());
+    for (int I=d.first(); I<=d.last(); ++I)
+      pa[idx++] = e(I);
+    VecRestoreArray(v, &pa);
+  }
+  template <class T>
+  static
+  void copy(const Engine<1, T, Brick>& e, Vec v)
+  {
+    PetscScalar *pa;
+    VecGetArray(v, &pa);
+    int idx=0;
+    Interval<1> d(e.domain());
+    for (int I=d.first(); I<=d.last(); ++I)
+      e(I) = pa[idx++];
+    VecRestoreArray(v, &pa);
+  }
+};
+
+template <>
+struct PoomaDACopy<2> {
+  template <class T>
+  static
+  void copy(Vec v, const Engine<2, T, Brick>& e)
+  {
+    PetscScalar *pa;
+    VecGetArray(v, &pa);
+    int idx=0;
+    Interval<2> d(e.domain());
+    for (int J=d[1].first(); J<=d[1].last(); ++J)
+      for (int I=d[0].first(); I<=d[0].last(); ++I)
+       pa[idx++] = e(I, J);
+    VecRestoreArray(v, &pa);
+  }
+  template <class T>
+  static
+  void copy(const Engine<2, T, Brick>& e, Vec v)
+  {
+    PetscScalar *pa;
+    VecGetArray(v, &pa);
+    int idx=0;
+    Interval<2> d(e.domain());
+    for (int J=d[1].first(); J<=d[1].last(); ++J)
+      for (int I=d[0].first(); I<=d[0].last(); ++I)
+       e(I, J) = pa[idx++];
+    VecRestoreArray(v, &pa);
+  }
+};
+
+template <>
+struct PoomaDACopy<3> {
+  template <class T>
+  static
+  void copy(Vec v, const Engine<3, T, Brick>& e)
+  {
+    PetscScalar *pa;
+    VecGetArray(v, &pa);
+    int idx=0;
+    Interval<3> d(e.domain());
+    for (int K=d[2].first(); K<=d[2].last(); ++K)
+      for (int J=d[1].first(); J<=d[1].last(); ++J)
+       for (int I=d[0].first(); I<=d[0].last(); ++I)
+         pa[idx++] = e(I, J, K);
+    VecRestoreArray(v, &pa);
+  }
+  template <class T>
+  static
+  void copy(const Engine<3, T, Brick>& e, Vec v)
+  {
+    PetscScalar *pa;
+    VecGetArray(v, &pa);
+    int idx=0;
+    Interval<3> d(e.domain());
+    for (int K=d[2].first(); K<=d[2].last(); ++K)
+      for (int J=d[1].first(); J<=d[1].last(); ++J)
+       for (int I=d[0].first(); I<=d[0].last(); ++I)
+         e(I, J, K) = pa[idx++];
+    VecRestoreArray(v, &pa);
+  }
+};
+
+
+
+/**
+ * Struct to wrap extra global information about a DA.
+ */
+
+struct PoomaDA {
+
+  /// Creates a PETSc DA from the specified layout.
+  /// Extra arguments are like DACreateNd, namely the periodicity
+  /// and stencil type and the stencil width.
+
+  template <class Layout>
+  PoomaDA(const Layout &l, DAPeriodicType pt, DAStencilType st, int sw);
+
+  ~PoomaDA()
+  {
+    delete[] info;
+    DADestroy(da);
+  }
+
+
+  /// Can use this as PETSc DA type.
+
+  operator DA() const { return da; }
+
+
+  /// Assign from POOMA engine to PETSc vector.
+
+  template <int Dim, class T, class EngineTag>
+  void assign(Vec v, const Engine<Dim, T, EngineTag> &e);
+
+  /// Assign from POOMA array to PETSc vector.
+
+  template <int Dim, class T, class EngineTag>
+  void assign(Vec v, const Array<Dim, T, EngineTag> &a)
+  {
+    this->assign(v, a.engine());
+  }
+
+  /// Assign from POOMA field to PETSc vector.
+
+  template <class MeshTag, class T, class EngineTag>
+  void assign(Vec v, const Field<MeshTag, T, EngineTag> &f)
+  {
+    this->assign(v, f.fieldEngine().engine());
+  }
+
+
+  /// Assign from PETSc vector to POOMA engine.
+
+  template <int Dim, class T, class EngineTag>
+  void assign(const Engine<Dim, T, EngineTag> &e, Vec v);
+
+  /// Assign from PETSc vector to POOMA array.
+
+  template <int Dim, class T, class EngineTag>
+  void assign(const Array<Dim, T, EngineTag> &a, Vec v)
+  {
+    this->assign(a.engine(), v);
+  }
+
+  /// Assign from PETSc vector to POOMA field.
+
+  template <class MeshTag, class T, class EngineTag>
+  void assign(const Field<MeshTag, T, EngineTag> &f, Vec v)
+  {
+    this->assign(f.fieldEngine().engine(), v);
+  }
+
+
+private:
+  DA da;
+  int n;
+  DALocalInfo *info;
+
+};
+
+
+template <class Layout>
+PoomaDA::PoomaDA(const Layout &l, DAPeriodicType pt, DAStencilType st, int sw)
+{
+  enum { Dim = Layout::dimensions };
+
+  // create DA
+  if (Dim == 1) {
+    DACreate1d(PETSC_COMM_WORLD,          /* MPI communicator */
+              pt,                        /* grid periodicity */
+              l.innerDomain()[0].size(), /* global domain size */
+              1,                         /* degrees of freedom */
+              sw,                        /* stencil width */
+              PETSC_NULL,                /* local domain sizes per dimension */
+              &this->da);
+  } else if (Dim == 2) {
+    DACreate2d(PETSC_COMM_WORLD,          /* MPI communicator */
+              pt,                        /* grid periodicity */
+              st,                        /* stencil type */
+              l.innerDomain()[0].size(),
+              l.innerDomain()[1].size(), /* global domain size */
+              PETSC_DECIDE, PETSC_DECIDE,/* # processors */
+              1,                         /* degrees of freedom */
+              sw,                        /* stencil width */
+              PETSC_NULL, PETSC_NULL,    /* local domain sizes per dimension */
+              &this->da);
+  } else if (Dim == 3) {
+    DACreate3d(PETSC_COMM_WORLD,          /* MPI communicator */
+              pt,                        /* grid periodicity */
+              st,                        /* stencil type */
+              l.innerDomain()[0].size(), l.innerDomain()[1].size(),
+              l.innerDomain()[2].size(), /* global domain size */
+              PETSC_DECIDE, PETSC_DECIDE,
+              PETSC_DECIDE,              /* # processors */
+              1,                         /* degrees of freedom */
+              sw,                        /* stencil width */
+              PETSC_NULL, PETSC_NULL,
+              PETSC_NULL,                /* local domain sizes per dimension */
+              &this->da);
+  }
+
+  // collect local information
+  int m, n, p;
+  DAGetInfo(this->da, PETSC_NULL, PETSC_NULL, PETSC_NULL, PETSC_NULL,
+            &m, &n, &p,
+            PETSC_NULL, PETSC_NULL, PETSC_NULL, PETSC_NULL);
+  this->n = m*n*p;
+  PInsist(Pooma::contexts() == this->n, "nr patches");
+  this->info = new DALocalInfo[this->n];
+  DAGetLocalInfo(this->da, &this->info[Pooma::context()]);
+
+  // distribute local information
+  // fixme - MPI_Allgather wrapper missing
+#if POOMA_MPI
+  MPI_Allgather(&this->info[Pooma::context()], sizeof(DALocalInfo), MPI_CHAR,
+               this->info, sizeof(DALocalInfo), MPI_CHAR,
+               MPI_COMM_WORLD);
+#endif
+}
+
+template <int Dim, class T, class EngineTag>
+void PoomaDA::assign(Vec v, const Engine<Dim, T, EngineTag> &e)
+{
+  typedef Engine<Dim, T, EngineTag> Engine_t;
+  typedef typename NewEngine<Engine_t, Interval<Dim> >::Type_t ViewEngine_t;
+
+  // our local brick engine
+  Engine<Dim, T, Brick> local_e;
+  Interval<Dim> local_I;
+
+  // loop over all DA patches
+  for (int i=0; i<this->n; ++i) {
+       // create DA patch context local pooma array
+        Interval<Dim> 
lPatch(PoomaDAGetDomain<Dim>::innerDomain(this->info[i]));
+       Array<Dim, T, Remote<Brick> > a;
+       a.engine() = Engine<Dim, T, Remote<Brick> >(i, lPatch);
+       Array<Dim, T, typename ViewEngine_t::Tag_t> e_array(ViewEngine_t(e, 
lPatch));
+       a = e_array;
+
+       // remember local engine
+       if (i == Pooma::context()) {
+         local_e = a.engine().localEngine();
+         local_I = lPatch;
+       }
+  }
+  Pooma::blockAndEvaluate();
+
+  // do the copy locally
+  PoomaDACopy<Dim>::copy(v, local_e);
+}
+
+template <int Dim, class T, class EngineTag>
+void PoomaDA::assign(const Engine<Dim, T, EngineTag> &e, Vec v)
+{
+  typedef Engine<Dim, T, EngineTag> Engine_t;
+  typedef typename NewEngine<Engine_t, Interval<Dim> >::Type_t ViewEngine_t;
+
+  // our local brick engine
+  Interval<Dim> 
local_I(PoomaDAGetDomain<Dim>::innerDomain(this->info[Pooma::context()]));
+  Engine<Dim, T, Brick> local_e(local_I);
+
+  // copy into the local brick
+  // if it were not the different array index ordering we could construct
+  // the brick engine with external data and avoid the double copying
+  PoomaDACopy<Dim>::copy(local_e, v);
+
+  // loop over all DA patches
+  for (int i=0; i<this->n; ++i) {
+       // create DA patch context local pooma array
+       Interval<Dim> lPatch(PoomaDAGetDomain<Dim>::innerDomain(this->info[i]));
+       Array<Dim, T, Remote<Brick> > a;
+       a.engine() = Engine<Dim, T, Remote<Brick> >(i, lPatch);
+
+       // override local engine with our one
+       if (Pooma::context() == i)
+         a.engine().localEngine() = local_e;
+
+       // distribute the copy
+       Array<Dim, T, typename ViewEngine_t::Tag_t> e_array;
+       e_array.engine() = ViewEngine_t(e, lPatch);
+       e_array = a;
+  }
+}
+
+
+} // namespace Pooma
+
+#endif

reply via email to

[Prev in Thread] Current Thread [Next in Thread]