Home Research Papers Posters Manuals Talks Download People Help Wiki Internal

libs/ck-libs/ParFUM/mesh.C

Go to the documentation of this file.

/*
Finite Element Method (FEM) Framework for Charm++
Parallel Programming Lab, Univ. of Illinois 2002
Orion Sky Lawlor, olawlor@acm.org, 12/20/2002

FEM Implementation file: mesh creation and user-data manipulation.
*/

#include "ParFUM.h"
#include "ParFUM_internals.h"

extern int femVersion;

/*** IDXL Interface ****/
FEM_Comm_Holder::FEM_Comm_Holder(FEM_Comm *sendComm, FEM_Comm *recvComm)
        :comm(sendComm,recvComm)
{
        registered=false;
        idx=-1;
}
void FEM_Comm_Holder::registerIdx(IDXL_Chunk *c) {
        assert(!registered);
        IDXL_Chunk *owner=c;
        if (idx!=-1) // had an old index: try to get it back
                idx=owner->addStatic(&comm,idx);
        else //No old index:
                idx=owner->addStatic(&comm);
}
void FEM_Comm_Holder::pup(PUP::er &p) {
        p|idx;
        if (p.isUnpacking() && idx!=-1)
        { // Try to grab the same index we had on our old processor:
                registerIdx(IDXL_Chunk::get("FEM_Comm_Holder::pup"));
        }
}

FEM_Comm_Holder::~FEM_Comm_Holder(void)
{
        if (registered)
        { // Try to unregister from IDXL:
                const char *caller="FEM_Comm_Holder::~FEM_Comm_Holder";
                IDXL_Chunk *owner=IDXL_Chunk::getNULL();
                if (owner) owner->destroy(idx,caller);
        }
}

/******* FEM_Mesh API ******/

static void checkIsSet(int fem_mesh,bool wantSet,const char *caller) {
        if (FEM_Mesh_is_set(fem_mesh)!=wantSet) {
                const char *msg=wantSet?"This mesh (%d) is not a setting mesh":
                        "This mesh (%d) is not a getting mesh";
                FEM_Abort(caller,msg,fem_mesh);
        }
}

/* Connectivity: Map calls to appropriate version of FEM_Mesh_data */
CDECL void
FEM_Mesh_conn(int fem_mesh,int entity,
        int *conn, int firstItem, int length, int width)
{
        FEM_Mesh_data(fem_mesh,entity,FEM_CONN, conn, firstItem,length, FEM_INDEX_0, width);
}
FDECL void
FTN_NAME(FEM_MESH_CONN,fem_mesh_conn)(int *fem_mesh,int *entity,
        int *conn, int *firstItem,int *length, int *width)
{
        //Can't just call the C version of this routine, because we use 1-based indices:
        FEM_Mesh_data(*fem_mesh,*entity,FEM_CONN, conn, *firstItem-1,*length, FEM_INDEX_1, *width);
}

CDECL void
FEM_Mesh_set_conn(int fem_mesh,int entity,
        const int *conn, int firstItem,int length, int width)
{
        checkIsSet(fem_mesh,true,"FEM_Mesh_set_conn");
        FEM_Mesh_conn(fem_mesh,entity,(int *)conn,firstItem,length,width);
}
FDECL void
FTN_NAME(FEM_MESH_SET_CONN,fem_mesh_set_conn)(int *fem_mesh,int *entity,
        const int *conn, int *firstItem,int *length, int *width)
{
        checkIsSet(*fem_mesh,true,"fem_mesh_set_conn");
        FTN_NAME(FEM_MESH_CONN,fem_mesh_conn)(fem_mesh,entity,(int *)conn,firstItem,length,width);
}

CDECL void
FEM_Mesh_get_conn(int fem_mesh,int entity,
        int *conn, int firstItem,int length, int width)
{
        checkIsSet(fem_mesh,false,"FEM_Mesh_get_conn");
        FEM_Mesh_conn(fem_mesh,entity,conn,firstItem,length,width);
}
FDECL void
FTN_NAME(FEM_MESH_GET_CONN,fem_mesh_get_conn)(int *fem_mesh,int *entity,
        int *conn, int *firstItem,int *length, int *width)
{
        checkIsSet(*fem_mesh,false,"fem_mesh_get_conn");
        FTN_NAME(FEM_MESH_CONN,fem_mesh_conn)(fem_mesh,entity,conn,firstItem,length,width);
}


/* Data: map to FEM_Mesh_offset */
CDECL void
FEM_Mesh_data(int fem_mesh,int entity,int attr,
        void *data, int firstItem,int length, int datatype,int width)
{
        IDXL_Layout lo(datatype,width);
        FEM_Mesh_data_layout(fem_mesh,entity,attr,data,firstItem,length,lo);
}

FORTRAN_AS_C(FEM_MESH_DATA,FEM_Mesh_data,fem_mesh_data,
        (int *fem_mesh,int *entity,int *attr,void *data,int *firstItem,int *length,int *datatype,int *width),
        (*fem_mesh,*entity,*attr,data,*firstItem-1,*length,*datatype,*width)
)


CDECL void
FEM_Mesh_set_data(int fem_mesh,int entity,int attr,
        const void *data, int firstItem,int length, int datatype,int width)
{
        checkIsSet(fem_mesh,true,"FEM_Mesh_set_data");
        FEM_Mesh_data(fem_mesh,entity,attr,(void *)data,firstItem,length,datatype,width);
}
FORTRAN_AS_C(FEM_MESH_SET_DATA,FEM_Mesh_set_data,fem_mesh_set_data,
        (int *fem_mesh,int *entity,int *attr,void *data,int *firstItem,int *length,int *datatype,int *width),
        (*fem_mesh,*entity,*attr,data,*firstItem-1,*length,*datatype,*width)
)

CDECL void
FEM_Mesh_get_data(int fem_mesh,int entity,int attr,
        void *data, int firstItem,int length, int datatype,int width)
{
        checkIsSet(fem_mesh,false,"FEM_Mesh_get_data");
        FEM_Mesh_data(fem_mesh,entity,attr,data,firstItem,length,datatype,width);
}
FORTRAN_AS_C(FEM_MESH_GET_DATA,FEM_Mesh_get_data,fem_mesh_get_data,
        (int *fem_mesh,int *entity,int *attr,void *data,int *firstItem,int *length,int *datatype,int *width),
        (*fem_mesh,*entity,*attr,data,*firstItem-1,*length,*datatype,*width)
)

CDECL void
FEM_Mesh_data_layout(int fem_mesh,int entity,int attr,
        void *data, int firstItem,int length, IDXL_Layout_t layout)
{
        const char *caller="FEM_Mesh_data_layout";
        FEM_Mesh_data_layout(fem_mesh,entity,attr,data,firstItem,length,
                IDXL_Layout_List::get().get(layout,caller));
}

FORTRAN_AS_C(FEM_MESH_DATA_LAYOUT,FEM_Mesh_data_layout,fem_mesh_data_layout,
        (int *fem_mesh,int *entity,int *attr,void *data,int *firstItem,int *length,int *layout),
        (*fem_mesh,*entity,*attr,data,*firstItem-1,*length,*layout)
)

CDECL void
FEM_Mesh_data_offset(int fem_mesh,int entity,int attr,
        void *data, int firstItem,int length,
        int type,int width, int offsetBytes,int distanceBytes,int skewBytes)
{
        const char *caller="FEM_Mesh_data_offset";
        FEM_Mesh_data_layout(fem_mesh,entity,attr,data,firstItem,length,
                IDXL_Layout(type,width,offsetBytes,distanceBytes,skewBytes));
}
FORTRAN_AS_C(FEM_MESH_DATA_OFFSET,FEM_Mesh_data_offset,fem_mesh_data_offset,
        (int *fem_mesh,int *entity,int *attr,
         void *data,int *firstItem,int *length,
         int *type,int *width,int *offset,int *distance,int *skew),
        (*fem_mesh,*entity,*attr,
         data,*firstItem-1,*length,
         *type,*width,*offset,*distance,*skew)
)


void FEM_Register_array(int fem_mesh,int entity,int attr,
        void *data, int datatype,int width,int firstItem){
        IDXL_Layout lo(datatype,width);
/*      if(attr == FEM_CONN){
                printf("CONN width %d \n",width);
                int len = FEM_Mesh_get_length(fem_mesh,entity);
                int *connd = (int *)data;
                for(int i=0;i<len;i++){
                        printf("%d -> (%d %d %d) \n",i+1,connd[3*i],connd[3*i+1],connd[3*i+2]);
                }
        }
        printf("firstItem %d \n",firstItem);*/
        FEM_Register_array_layout(fem_mesh,entity,attr,data,firstItem,lo);
}

void FEM_Register_array_layout(int fem_mesh,int entity,int attr,
        void *data, IDXL_Layout_t layout,int firstItem){
        const char *caller="FEM_Register_array_layout";
        FEM_Register_array_layout(fem_mesh,entity,attr,data,firstItem,
                IDXL_Layout_List::get().get(layout,caller));

}

/*registration api */
CDECL void
FEM_Register_array(int fem_mesh,int entity,int attr,
        void *data, int datatype,int width)
{
        FEM_Register_array(fem_mesh,entity,attr,data,datatype,width,0);
}

CDECL void
FEM_Register_array_layout(int fem_mesh,int entity,int attr,
        void *data, IDXL_Layout_t layout){
        FEM_Register_array_layout(fem_mesh,entity,attr,data,layout,0);
}


CDECL void
FEM_Register_entity(int fem_mesh,int entity,void *data,
                int len,int max,FEM_Mesh_alloc_fn fn) {
                FEM_Register_entity_impl(fem_mesh,entity,data,len,max,fn);
}

FORTRAN_AS_C(FEM_REGISTER_ARRAY,FEM_Register_array,fem_register_array,
        (int *fem_mesh,int *entity,int *attr,void *data,int *datatype,int *width),(*fem_mesh,*entity,*attr,data,*datatype,*width,0))


FORTRAN_AS_C(FEM_REGISTER_ARRAY_LAYOUT,FEM_Register_array_layout,fem_register_array_layout,
        (int *fem_mesh,int *entity,int *attr,void *data,int *layout),(*fem_mesh,*entity,*attr,data,*layout,0))

FORTRAN_AS_C(FEM_REGISTER_ENTITY,FEM_Register_entity,fem_register_entity,
        (int *fem_mesh,int *entity,void *data,int *len,int *max,FEM_Mesh_alloc_fn fn),(*fem_mesh,*entity,data,*len,*max,fn))


// User data API:
CDECL void
FEM_Mesh_pup(int fem_mesh,int dataTag,FEM_Userdata_fn fn,void *data) {
        const char *caller="FEM_Mesh_pup"; FEMAPI(caller);
        FEM_Mesh *m=FEM_Mesh_lookup(fem_mesh,caller);
        FEM_Userdata_item &i=m->udata.find(dataTag);
        FEM_Userdata_pupfn f(fn,data);
        if (m->isSetting()) i.store(f);
        else /* m->isGetting() */ {
                if (!i.hasStored())
                        FEM_Abort(caller,"Never stored any user data at tag %d",dataTag);
                i.restore(f);
        }
}
FORTRAN_AS_C(FEM_MESH_PUP,FEM_Mesh_pup,fem_mesh_pup,
        (int *m,int *t,FEM_Userdata_fn fn,void *data), (*m,*t,fn,data))

// Accessor API:
CDECL void
FEM_Mesh_set_length(int fem_mesh,int entity,int newLength) {
        const char *caller="FEM_Mesh_set_length"; FEMAPI(caller);
        checkIsSet(fem_mesh,true,caller);
        FEM_Entity_lookup(fem_mesh,entity,caller)->setLength(newLength);
}
FORTRAN_AS_C(FEM_MESH_SET_LENGTH,FEM_Mesh_set_length,fem_mesh_set_length,
        (int *fem_mesh,int *entity,int *newLength),
        (*fem_mesh,*entity,*newLength)
)


CDECL int
FEM_Mesh_get_length(int fem_mesh,int entity) {
        const char *caller="FEM_Mesh_get_length"; FEMAPI(caller);
        int len=FEM_Entity_lookup(fem_mesh,entity,caller)->size();
        return len;
}
FORTRAN_AS_C_RETURN(int,
        FEM_MESH_GET_LENGTH,FEM_Mesh_get_length,fem_mesh_get_length,
        (int *fem_mesh,int *entity),(*fem_mesh,*entity)
)


CDECL void
FEM_Mesh_set_width(int fem_mesh,int entity,int attr,int newWidth) {
        const char *caller="FEM_Mesh_set_width";
        FEMAPI(caller);
        checkIsSet(fem_mesh,true,caller);
        FEM_Attribute_lookup(fem_mesh,entity,attr,caller)->setWidth(newWidth,caller);
}
FORTRAN_AS_C(FEM_MESH_SET_WIDTH,FEM_Mesh_set_width,fem_mesh_set_width,
        (int *fem_mesh,int *entity,int *attr,int *newWidth),
        (*fem_mesh,*entity,*attr,*newWidth)
)

CDECL int
FEM_Mesh_get_width(int fem_mesh,int entity,int attr) {
        const char *caller="FEM_Mesh_get_width";
        FEMAPI(caller);
        return FEM_Attribute_lookup(fem_mesh,entity,attr,caller)->getWidth();
}
FORTRAN_AS_C_RETURN(int,
        FEM_MESH_GET_WIDTH,FEM_Mesh_get_width,fem_mesh_get_width,
        (int *fem_mesh,int *entity,int *attr),(*fem_mesh,*entity,*attr)
)

CDECL int
FEM_Mesh_get_datatype(int fem_mesh,int entity,int attr) {
        const char *caller="FEM_Mesh_get_datatype";
        FEMAPI(caller);
        return FEM_Attribute_lookup(fem_mesh,entity,attr,caller)->getDatatype();
}
FORTRAN_AS_C_RETURN(int,
        FEM_MESH_GET_DATATYPE,FEM_Mesh_get_datatype,fem_mesh_get_datatype,
        (int *fem_mesh,int *entity,int *attr),(*fem_mesh,*entity,*attr)
)

CDECL int
FEM_Mesh_is_set(int fem_mesh) /* return 1 if this is a writing mesh */
{
        return (FEM_Mesh_lookup(fem_mesh,"FEM_Mesh_is_get")->isSetting())?1:0;
}
FORTRAN_AS_C_RETURN(int,
        FEM_MESH_IS_SET,FEM_Mesh_is_set,fem_mesh_is_set,
        (int *fem_mesh),(*fem_mesh)
)

CDECL int
FEM_Mesh_is_get(int fem_mesh) /* return 1 if this is a readable mesh */
{
        return (!FEM_Mesh_lookup(fem_mesh,"FEM_Mesh_is_get")->isSetting())?1:0;
}
FORTRAN_AS_C_RETURN(int,
        FEM_MESH_IS_GET,FEM_Mesh_is_get,fem_mesh_is_get,
        (int *fem_mesh),(*fem_mesh)
)

CDECL void
FEM_Mesh_become_get(int fem_mesh) /* Make this a readable mesh */
{ FEM_Mesh_lookup(fem_mesh,"FEM_Mesh_become_get")->becomeGetting(); }
FORTRAN_AS_C(FEM_MESH_BECOME_GET,FEM_Mesh_become_get,fem_mesh_become_get, (int *m),(*m))

CDECL void
FEM_Mesh_become_set(int fem_mesh)
{ FEM_Mesh_lookup(fem_mesh,"FEM_Mesh_become_get")->becomeSetting(); }
FORTRAN_AS_C(FEM_MESH_BECOME_SET,FEM_Mesh_become_set,fem_mesh_become_set, (int *m),(*m))


CDECL IDXL_t
FEM_Comm_shared(int fem_mesh,int entity) {
        const char *caller="FEM_Comm_shared";
        FEMAPI(caller);
        if (entity!=FEM_NODE) FEM_Abort(caller,"Only shared nodes supported");
        return FEM_Mesh_lookup(fem_mesh,caller)->node.
                sharedIDXL.getIndex(IDXL_Chunk::get(caller));
}
FORTRAN_AS_C_RETURN(int,
        FEM_COMM_SHARED,FEM_Comm_shared,fem_comm_shared,
        (int *fem_mesh,int *entity),(*fem_mesh,*entity)
)

CDECL IDXL_t
FEM_Comm_ghost(int fem_mesh,int entity) {
        const char *caller="FEM_Comm_ghost";
        FEMAPI(caller);
        FEM_Entity *e=FEM_Entity_lookup(fem_mesh,entity,caller);
        if (e->isGhost()) FEM_Abort(caller,"Can only call FEM_Comm_ghost on real entity type");
        return e->ghostIDXL.getIndex(IDXL_Chunk::get(caller));
}
FORTRAN_AS_C_RETURN(int,
        FEM_COMM_GHOST,FEM_Comm_ghost,fem_comm_ghost,
        (int *fem_mesh,int *entity),(*fem_mesh,*entity)
)


// Internal API:
void FEM_Mesh_data_layout(int fem_mesh,int entity,int attr,
        void *data, int firstItem,int length, const IDXL_Layout &layout)
{
        if (femVersion == 0 && length==0) return;
        const char *caller="FEM_Mesh_data";
        FEMAPI(caller);
        FEM_Mesh *m=FEM_Mesh_lookup(fem_mesh,caller);
        FEM_Attribute *a=m->lookup(entity,caller)->lookup(attr,caller);

        if (m->isSetting())
                a->set(data,firstItem,length,layout,caller);
        else /* m->isGetting()*/
                a->get(data,firstItem,length,layout,caller);
}

void FEM_Register_array_layout(int fem_mesh,int entity,int attr,void *data,int firstItem,const IDXL_Layout &layout){
        const char *caller="FEM_Register_array";
        FEMAPI(caller);
        FEM_Mesh *m=FEM_Mesh_lookup(fem_mesh,caller);
        FEM_Entity *e = m->lookup(entity,caller);
        int length = e->size();
        //should actually be a call on the entity
        int max = e->getMax();
        FEM_Attribute *a = e->lookup(attr,caller);
        if(a->getWidth() == 0){
        //      a->setDatatype(layout.type);
                a->setWidth(layout.width);
        }


        if(m->isSetting()){
          a->register_data(data,length,max,layout,caller);
        }else{
          a->get(data,firstItem,length,layout,caller);
          //replace the attribute's data array with the user's data
          a->register_data(data,max/*length*/,max,layout,caller);
        }
}
void FEM_Register_entity_impl(int fem_mesh,int entity,void *args,int len,int max,FEM_Mesh_alloc_fn fn){
        char *caller = "FEM_Register_entity";
        FEM_Mesh *m=FEM_Mesh_lookup(fem_mesh,caller);
/*      if(!m->isSetting()){
                CmiAbort("Register entity called on mesh that can't be written into");
        }
*/
        FEM_Entity *e = m->lookup(entity,caller);
        e->setMaxLength(len,max,args,fn);
}

FEM_Entity *FEM_Entity_lookup(int fem_mesh,int entity,const char *caller) {
        return FEM_Mesh_lookup(fem_mesh,caller)->lookup(entity,caller);
}
FEM_Attribute *FEM_Attribute_lookup(int fem_mesh,int entity,int attr,const char *caller) {
        return FEM_Entity_lookup(fem_mesh,entity,caller)->lookup(attr,caller);
}

CDECL int FEM_Mesh_get_entities(int fem_mesh, int *entities) {
        const char *caller="FEM_Mesh_get_entities";
        FEMAPI(caller);
        return FEM_Mesh_lookup(fem_mesh,caller)->getEntities(entities);
}
FORTRAN_AS_C_RETURN(int,
        FEM_MESH_GET_ENTITIES,FEM_Mesh_get_entities,fem_mesh_get_entities,
        (int *mesh,int *ent), (*mesh,ent)
)

CDECL int FEM_Mesh_get_attributes(int fem_mesh,int entity,int *attributes) {
        const char *caller="FEM_Mesh_get_attributes";
        FEMAPI(caller);
        return FEM_Entity_lookup(fem_mesh,entity,caller)->getAttrs(attributes);
}
FORTRAN_AS_C_RETURN(int,
        FEM_MESH_GET_ATTRIBUTES,FEM_Mesh_get_attributes,fem_mesh_get_attributes,
        (int *mesh,int *ent,int *attrs), (*mesh,*ent,attrs)
)

/************** FEM_Attribute ****************/

CDECL const char *FEM_Get_datatype_name(int datatype,char *storage) {
        switch(datatype) {
        case FEM_BYTE: return "FEM_BYTE";
        case FEM_INT: return "FEM_INT";
        case FEM_FLOAT: return "FEM_FLOAT";
        case FEM_DOUBLE: return "FEM_DOUBLE";
        case FEM_INDEX_0: return "FEM_INDEX_0";
        case FEM_INDEX_1: return "FEM_INDEX_1";
        };
        sprintf(storage,"unknown datatype code (%d)",datatype);
        return storage;
}

CDECL const char *FEM_Get_attr_name(int attr,char *storage)
{
        if (attr<FEM_ATTRIB_TAG_MAX)
        { //It's a user tag:
                sprintf(storage,"FEM_DATA+%d",attr-FEM_DATA);
                return storage;
        }
        switch(attr) {
        case FEM_CONN: return "FEM_CONN";
        case FEM_SPARSE_ELEM: return "FEM_SPARSE_ELEM";
        case FEM_COOR: return "FEM_COOR";
        case FEM_GLOBALNO: return "FEM_GLOBALNO";
        case FEM_PARTITION: return "FEM_PARTITION";
        case FEM_SYMMETRIES: return "FEM_SYMMETRIES";
        case FEM_NODE_PRIMARY: return "FEM_NODE_PRIMARY";
        case FEM_NODE_ELEM_ADJACENCY: return "FEM_NODE_ELEM_ADJACENCY";
        case FEM_NODE_NODE_ADJACENCY: return "FEM_NODE_NODE_ADJACENCY";
        case FEM_ELEM_ELEM_ADJACENCY: return "FEM_ELEM_ELEM_ADJACENCY";
        case FEM_ELEM_ELEM_ADJ_TYPES: return "FEM_ELEM_ELEM_ADJ_TYPES";
        case FEM_IS_VALID_ATTR: return "FEM_IS_VALID_ATTR";
        case FEM_MESH_SIZING: return "FEM_MESH_SIZING";

        default: break;
        };
        sprintf(storage,"unknown attribute code (%d)",attr);
        return storage;
}

//Abort with a nice error message saying:
// Our <field> was previously set to <cur>; it cannot now be <operation> <next>
void FEM_Attribute::bad(const char *field,bool forRead,int cur,int next,const char *caller) const
{
        char nameStorage[256];
        const char *name=FEM_Get_attr_name(attr,nameStorage);
        char errBuf[1024];
        const char *cannotBe=NULL;
        if (forRead) {
                if (cur==-1) {
                        sprintf(errBuf,"The %s %s %s was never set-- it cannot now be read",
                                e->getName(),name,field);
                }
                else /* already had a value */
                        cannotBe="read as";
        }
        else /* for write */ {
                cannotBe="set to";
        }
        if (cannotBe!=NULL) /* Use standard ... <something> cannot be <something>'d... error message */
                sprintf(errBuf,"The %s %s %s was previously set to %d; it cannot now be %s %d",
                        e->getName(),name,field,cur,cannotBe,next);

        FEM_Abort(caller,errBuf);
}


FEM_Attribute::FEM_Attribute(FEM_Entity *e_,int attr_)
                :e(e_),ghost(0),attr(attr_),width(0),datatype(-1), allocated(false)
{
        tryAllocate();
        if (femVersion == 0) width=-1;
}
void FEM_Attribute::pup(PUP::er &p) {
        // e, attr, and ghost are always set by the constructor
        p|width;
        if (p.isUnpacking() && femVersion > 0 && width<0)  width=0;
        p|datatype;
        if (p.isUnpacking()) tryAllocate();
}
void FEM_Attribute::pupSingle(PUP::er &p, int pupindx) {
  return;
}
FEM_Attribute::~FEM_Attribute() {}

void FEM_Attribute::setLength(int next,const char *caller) {
        int cur=getLength();
        if (next==cur) return; //Already set--nothing to do
        if (cur>0) bad("length",false,cur,next, caller);
        e->setLength(next);
        tryAllocate();
}

void FEM_Attribute::setWidth(int next,const char *caller) {
        int cur=getWidth();
        if (next==cur) return; //Already set--nothing to do
        if (cur>0) bad("width",false,cur,next, caller);
        width=next;
        tryAllocate();
        if (ghost) ghost->setWidth(width,caller);
}

void FEM_Attribute::setDatatype(int next,const char *caller) {
        int cur=getDatatype();
        if (next==cur) return; //Already set--nothing to do
        if (cur!=-1) bad("datatype",false,cur,next, caller);
        datatype=next;
        tryAllocate();
        if (ghost) ghost->setDatatype(datatype,caller);
}

void FEM_Attribute::copyShape(const FEM_Attribute &src) {
        setWidth(src.getWidth());
        if (src.getDatatype()!=-1)
          setDatatype(src.getDatatype()); //Automatically calls tryAllocate
}
void FEM_Attribute::set(const void *src, int firstItem,int length,
                const IDXL_Layout &layout, const char *caller)
{
        if (firstItem!=0) { /* If this isn't the start... */
                if (length!=1) /* And we're not setting one at a time */
                        CmiAbort("FEM_Mesh_data: unexpected firstItem");
        }

        if (femVersion == 0 && getRealLength() == -1) setLength(length);
        else if (getLength()==0) setLength(length);
        else if (length!=1 && length!=getLength())
                bad("length",false,getLength(),length, caller);

        int width=layout.width;
        if (femVersion==0 && getRealWidth()==-1) setWidth(width);
        else if (getWidth()==0) setWidth(width);
        else if (width!=getWidth())
                bad("width",false,getWidth(),width, caller);

        int datatype=layout.type;
        if (getDatatype()==-1) setDatatype(datatype);
        else if (datatype!=getDatatype())
                bad("datatype",false,getDatatype(),datatype, caller);

        /* Assert: our storage should be allocated now.
           Our subclass will actually copy user data */
}

void FEM_Attribute::get(void *dest, int firstItem,int length,
                const IDXL_Layout &layout, const char *caller)  const
{
        if (length==0) return; //Nothing to get
        if (length!=1 && length!=getLength())
                bad("length",true,getLength(),length, caller);

        int width=layout.width;
        if (width!=getWidth())
                bad("width",true,getWidth(),width, caller);

        int datatype=layout.type;
        if (datatype!=getDatatype())
                bad("datatype",true,getDatatype(),datatype, caller);

        /* our subclass will actually copy into user data */
}

/*check if the layout is the same as earlier */

void FEM_Attribute::register_data(void *user, int length,int max,
        const IDXL_Layout &layout, const char *caller){

                int width=layout.width;
                if (femVersion == 0 && getRealWidth()==-1) setWidth(width);
                else if (getWidth()==0){
                        setWidth(width);
                }else{
                        if (width!=getWidth()){
                                bad("width",false,getWidth(),width, caller);
                        }
                }

                int datatype=layout.type;
                if (getDatatype()==-1){
                        setDatatype(datatype);
                }else{
                        if (datatype!=getDatatype()){
                                bad("datatype",false,getDatatype(),datatype, caller);
                        }
                }

}

//Check if all three of length, width, and datatype are set.
// If so, call allocate.
void FEM_Attribute::tryAllocate(void) {
        int lenNull, widthNull;
        if (femVersion == 0) {
          // version 0 takes -1 as empty
          lenNull = (getRealLength()==-1);
          widthNull = (getRealWidth()==-1);
        }
        else {
          lenNull = (getLength()==0);
          widthNull = (getWidth()==0);
        }
        if ((!allocated) && !lenNull && !widthNull && getDatatype()!=-1) {
          allocated=true;
                allocate(getMax(),getWidth(),getDatatype());
        }
}

/*********************** DataAttribute *******************/
FEM_DataAttribute::FEM_DataAttribute(FEM_Entity *e,int myAttr)
        :FEM_Attribute(e,myAttr),
         char_data(0),int_data(0),float_data(0),double_data(0)
{
}
void FEM_DataAttribute::pup(PUP::er &p) {
        super::pup(p);
        switch(getDatatype()) {
        case -1: /* not allocated yet */ break;
        case FEM_BYTE:
          if (char_data) {
            /*if(p.isSizing()) {
              char_data->setRowLen(getEntity()->size());
              }*/
            char_data->pup(p);
          }
          break;
        case FEM_INT:
          if (int_data) {
            /*if(p.isSizing()) {
              int_data->setRowLen(getEntity()->size());
              }*/
            int_data->pup(p);
          }
          break;
        case FEM_FLOAT:
          if (float_data) {
            /*if(p.isSizing()) {
              float_data->setRowLen(getEntity()->size());
              }*/
            float_data->pup(p);
          }
          break;
        case FEM_DOUBLE:
          if (double_data) {
            /*if(p.isSizing()) {
              double_data->setRowLen(getEntity()->size());
              }*/
            double_data->pup(p);
          }
          break;
        default: CkAbort("Invalid datatype in FEM_DataAttribute::pup");
        }
}
void FEM_DataAttribute::pupSingle(PUP::er &p, int pupindx) {
        super::pupSingle(p,pupindx);
        switch(getDatatype()) {
        case -1: /* not allocated yet */ break;
        case FEM_BYTE:   if (char_data) char_data->pupSingle(p,pupindx); break;
        case FEM_INT:    if (int_data) int_data->pupSingle(p,pupindx); break;
        case FEM_FLOAT:  if (float_data) float_data->pupSingle(p,pupindx); break;
        case FEM_DOUBLE: if (double_data) double_data->pupSingle(p,pupindx); break;
        default: CkAbort("Invalid datatype in FEM_DataAttribute::pupSingle");
        }
}
FEM_DataAttribute::~FEM_DataAttribute() {
        if (char_data) delete char_data;
        if (int_data) delete int_data;
        if (float_data) delete float_data;
        if (double_data) delete double_data;

}

template <class T>
inline void setTableData(const void *user, int firstItem, int length,
        IDXL_LAYOUT_PARAM, AllocTable2d<T> *table)
{
        for (int r=0;r<length;r++) {
                register T *tableRow=table->getRow(firstItem+r);
                for (int c=0;c<width;c++)
                        tableRow[c]=IDXL_LAYOUT_DEREF(T,user,r,c);
        }
}

template <class T>
inline void getTableData(void *user, int firstItem, int length,
        IDXL_LAYOUT_PARAM, const AllocTable2d<T> *table)
{
        for (int r=0;r<length;r++) {
                register const T *tableRow=table->getRow(firstItem+r);
                for (int c=0;c<width;c++)
                        IDXL_LAYOUT_DEREF(T,user,r,c)=tableRow[c];
        }

}

void FEM_DataAttribute::set(const void *u, int f,int l,