/[escript]/trunk/escript/src/DataTagged.cpp
ViewVC logotype

Diff of /trunk/escript/src/DataTagged.cpp

Parent Directory Parent Directory | Revision Log Revision Log | View Patch Patch

revision 594 by gross, Tue Mar 14 23:36:08 2006 UTC revision 1628 by phornby, Fri Jul 11 13:12:46 2008 UTC
# Line 1  Line 1 
 // $Id$  
1    
2  /*  /* $Id$ */
3   ******************************************************************************  
4   *                                                                            *  /*******************************************************
5   *       COPYRIGHT  ACcESS 2004 -  All Rights Reserved                        *   *
6   *                                                                            *   *           Copyright 2003-2007 by ACceSS MNRF
7   * This software is the property of ACcESS. No part of this code              *   *       Copyright 2007 by University of Queensland
8   * may be copied in any form or by any means without the expressed written    *   *
9   * consent of ACcESS.  Copying, use or modification of this software          *   *                http://esscc.uq.edu.au
10   * by any unauthorised person is illegal unless that person has a software    *   *        Primary Business: Queensland, Australia
11   * license agreement with ACcESS.                                             *   *  Licensed under the Open Software License version 3.0
12   *                                                                            *   *     http://www.opensource.org/licenses/osl-3.0.php
13   ******************************************************************************   *
14  */   *******************************************************/
15    
16  #include "DataTagged.h"  #include "DataTagged.h"
17    #include "esysUtils/esys_malloc.h"
18    
19  #include "DataConstant.h"  #include "DataConstant.h"
20  #include "DataException.h"  #include "DataException.h"
21    #ifdef USE_NETCDF
22    #include <netcdfcpp.h>
23    #endif
24    
25  using namespace std;  using namespace std;
26    
# Line 80  DataTagged::DataTagged(const FunctionSpa Line 83  DataTagged::DataTagged(const FunctionSpa
83    }    }
84  }  }
85    
86    DataTagged::DataTagged(const FunctionSpace& what,
87                           const DataArrayView::ShapeType &shape,
88                           const TagListType& tags,
89                           const ValueType& data)
90      : DataAbstract(what)
91    {
92      // alternative constructor
93      // not unit_tested tested yet
94    
95      // copy the data
96      m_data=data;
97    
98      // create the view of the data
99      DataArrayView tempView(m_data,shape);
100      setPointDataView(tempView);
101    
102      // create the tag lookup map
103      for (int sampleNo=0; sampleNo<getNumSamples(); sampleNo++) {
104        m_offsetLookup.insert(DataMapType::value_type(sampleNo,tags[sampleNo]));
105      }
106    }
107    
108    
109  DataTagged::DataTagged(const DataTagged& other)  DataTagged::DataTagged(const DataTagged& other)
110    : DataAbstract(other.getFunctionSpace()),    : DataAbstract(other.getFunctionSpace()),
111    m_data(other.m_data),    m_data(other.m_data),
# Line 110  DataTagged::DataTagged(const DataConstan Line 136  DataTagged::DataTagged(const DataConstan
136    setPointDataView(temp);    setPointDataView(temp);
137  }  }
138    
 void  
 DataTagged::reshapeDataPoint(const DataArrayView::ShapeType& shape)  
 {  
   // can only reshape a rank zero data point  
   if (getPointDataView().getRank()!=0) {  
     stringstream temp;  
     temp << "Error - Can only reshape Data with data points of rank 0. "  
      << "This Data has data points with rank: "  
      << getPointDataView().getRank();  
     throw DataException(temp.str());  
   }  
   
   // allocate enough space for all values  
   DataArrayView::ValueType newData(DataArrayView::noValues(shape)*(m_offsetLookup.size()+1));  
   DataArrayView newView(newData,shape);  
   newView.copy(0,getDefaultValue()());  
   
   // loop through the tag values  
   DataMapType::iterator pos;  
   DataArrayView::ValueType::size_type tagOffset=DataArrayView::noValues(shape);  
   for (pos=m_offsetLookup.begin();pos!=m_offsetLookup.end();pos++){  
     newView.copy(tagOffset,m_data[pos->second]);  
     pos->second=tagOffset;  
     tagOffset+=DataArrayView::noValues(shape);  
   }  
   m_data=newData;  
   DataArrayView temp(m_data,shape);  
   setPointDataView(temp);  
 }  
   
139  DataAbstract*  DataAbstract*
140  DataTagged::getSlice(const DataArrayView::RegionType& region) const  DataTagged::getSlice(const DataArrayView::RegionType& region) const
141  {  {
# Line 200  DataTagged::setSlice(const DataAbstract* Line 196  DataTagged::setSlice(const DataAbstract*
196    if (getPointDataView().getRank()!=region.size()) {    if (getPointDataView().getRank()!=region.size()) {
197      throw DataException("Error - Invalid slice region.");      throw DataException("Error - Invalid slice region.");
198    }    }
199    if (otherTemp->getPointDataView().getRank()>0 and !other->getPointDataView().checkShape(regionShape)) {    if (otherTemp->getPointDataView().getRank()>0 && !other->getPointDataView().checkShape(regionShape)) {
200      throw DataException (other->getPointDataView().createShapeErrorMessage(      throw DataException (other->getPointDataView().createShapeErrorMessage(
201                           "Error - Couldn't copy slice due to shape mismatch.",regionShape));                           "Error - Couldn't copy slice due to shape mismatch.",regionShape));
202    }    }
# Line 305  DataTagged::addTaggedValues(const TagLis Line 301  DataTagged::addTaggedValues(const TagLis
301         << " doesn't match number of values: " << values.size();         << " doesn't match number of values: " << values.size();
302        throw DataException(temp.str());        throw DataException(temp.str());
303      } else {      } else {
304        for (int i=0;i<tagKeys.size();i++) {        unsigned int i;
305          for (i=0;i<tagKeys.size();i++) {
306          addTaggedValue(tagKeys[i],values[i]);          addTaggedValue(tagKeys[i],values[i]);
307        }        }
308      }      }
# Line 422  DataTagged::extractData(ifstream& archiv Line 419  DataTagged::extractData(ifstream& archiv
419    return(m_data.extractData(archiveFile, noValues));    return(m_data.extractData(archiveFile, noValues));
420  }  }
421  void  void
422    DataTagged::symmetric(DataAbstract* ev)
423    {
424      DataTagged* temp_ev=dynamic_cast<DataTagged*>(ev);
425      if (temp_ev==0) {
426        throw DataException("Error - DataTagged::symmetric casting to DataTagged failed (probably a programming error).");
427      }
428      const DataTagged::DataMapType& thisLookup=getTagLookup();
429      DataTagged::DataMapType::const_iterator i;
430      DataTagged::DataMapType::const_iterator thisLookupEnd=thisLookup.end();
431      for (i=thisLookup.begin();i!=thisLookupEnd;i++) {
432          temp_ev->addTaggedValue(i->first,temp_ev->getDefaultValue());
433          DataArrayView thisView=getDataPointByTag(i->first);
434          DataArrayView evView=temp_ev->getDataPointByTag(i->first);
435          DataArrayView::symmetric(thisView,0,evView,0);
436      }
437      DataArrayView::symmetric(getDefaultValue(),0,temp_ev->getDefaultValue(),0);
438    }
439    void
440    DataTagged::nonsymmetric(DataAbstract* ev)
441    {
442      DataTagged* temp_ev=dynamic_cast<DataTagged*>(ev);
443      if (temp_ev==0) {
444        throw DataException("Error - DataTagged::nonsymmetric casting to DataTagged failed (probably a programming error).");
445      }
446      const DataTagged::DataMapType& thisLookup=getTagLookup();
447      DataTagged::DataMapType::const_iterator i;
448      DataTagged::DataMapType::const_iterator thisLookupEnd=thisLookup.end();
449      for (i=thisLookup.begin();i!=thisLookupEnd;i++) {
450          temp_ev->addTaggedValue(i->first,temp_ev->getDefaultValue());
451          DataArrayView thisView=getDataPointByTag(i->first);
452          DataArrayView evView=temp_ev->getDataPointByTag(i->first);
453          DataArrayView::nonsymmetric(thisView,0,evView,0);
454      }
455      DataArrayView::nonsymmetric(getDefaultValue(),0,temp_ev->getDefaultValue(),0);
456    }
457    void
458    DataTagged::trace(DataAbstract* ev, int axis_offset)
459    {
460      DataTagged* temp_ev=dynamic_cast<DataTagged*>(ev);
461      if (temp_ev==0) {
462        throw DataException("Error - DataTagged::trace casting to DataTagged failed (probably a programming error).");
463      }
464      const DataTagged::DataMapType& thisLookup=getTagLookup();
465      DataTagged::DataMapType::const_iterator i;
466      DataTagged::DataMapType::const_iterator thisLookupEnd=thisLookup.end();
467      for (i=thisLookup.begin();i!=thisLookupEnd;i++) {
468          temp_ev->addTaggedValue(i->first,temp_ev->getDefaultValue());
469          DataArrayView thisView=getDataPointByTag(i->first);
470          DataArrayView evView=temp_ev->getDataPointByTag(i->first);
471          DataArrayView::trace(thisView,0,evView,0, axis_offset);
472      }
473      DataArrayView::trace(getDefaultValue(),0,temp_ev->getDefaultValue(),0,axis_offset);
474    }
475    
476    void
477    DataTagged::transpose(DataAbstract* ev, int axis_offset)
478    {
479      DataTagged* temp_ev=dynamic_cast<DataTagged*>(ev);
480      if (temp_ev==0) {
481        throw DataException("Error - DataTagged::transpose casting to DataTagged failed (probably a programming error).");
482      }
483      const DataTagged::DataMapType& thisLookup=getTagLookup();
484      DataTagged::DataMapType::const_iterator i;
485      DataTagged::DataMapType::const_iterator thisLookupEnd=thisLookup.end();
486      for (i=thisLookup.begin();i!=thisLookupEnd;i++) {
487          temp_ev->addTaggedValue(i->first,temp_ev->getDefaultValue());
488          DataArrayView thisView=getDataPointByTag(i->first);
489          DataArrayView evView=temp_ev->getDataPointByTag(i->first);
490          DataArrayView::transpose(thisView,0,evView,0, axis_offset);
491      }
492      DataArrayView::transpose(getDefaultValue(),0,temp_ev->getDefaultValue(),0,axis_offset);
493    }
494    
495    void
496    DataTagged::swapaxes(DataAbstract* ev, int axis0, int axis1)
497    {
498      DataTagged* temp_ev=dynamic_cast<DataTagged*>(ev);
499      if (temp_ev==0) {
500        throw DataException("Error - DataTagged::swapaxes casting to DataTagged failed (probably a programming error).");
501      }
502      const DataTagged::DataMapType& thisLookup=getTagLookup();
503      DataTagged::DataMapType::const_iterator i;
504      DataTagged::DataMapType::const_iterator thisLookupEnd=thisLookup.end();
505      for (i=thisLookup.begin();i!=thisLookupEnd;i++) {
506          temp_ev->addTaggedValue(i->first,temp_ev->getDefaultValue());
507          DataArrayView thisView=getDataPointByTag(i->first);
508          DataArrayView evView=temp_ev->getDataPointByTag(i->first);
509          DataArrayView::swapaxes(thisView,0,evView,0,axis0,axis1);
510      }
511      DataArrayView::swapaxes(getDefaultValue(),0,temp_ev->getDefaultValue(),0,axis0,axis1);
512    }
513    
514    void
515  DataTagged::eigenvalues(DataAbstract* ev)  DataTagged::eigenvalues(DataAbstract* ev)
516  {  {
517    DataTagged* temp_ev=dynamic_cast<DataTagged*>(ev);    DataTagged* temp_ev=dynamic_cast<DataTagged*>(ev);
# Line 469  DataTagged::eigenvalues_and_eigenvectors Line 559  DataTagged::eigenvalues_and_eigenvectors
559    
560  }  }
561    
562    void
563    DataTagged::setToZero(){
564        DataArrayView::ValueType::size_type n=m_data.size();
565        for (int i=0; i<n ;++i) m_data[i]=0.;
566    }
567    
568    void
569    DataTagged::dump(const std::string fileName) const
570    {
571       #ifdef PASO_MPI
572       throw DataException("Error - DataTagged:: dump is not implemented for MPI yet.");
573       #endif
574       #ifdef USE_NETCDF
575       const int ldims=DataArrayView::maxRank+1;
576       const NcDim* ncdims[ldims];
577       NcVar *var, *tags_var;
578       int rank = getPointDataView().getRank();
579       int type=  getFunctionSpace().getTypeCode();
580       int ndims =0;
581       long dims[ldims];
582       const double* d_ptr=&(m_data[0]);
583       DataArrayView::ShapeType shape = getPointDataView().getShape();
584    
585       // netCDF error handler
586       NcError err(NcError::verbose_nonfatal);
587       // Create the file.
588       NcFile dataFile(fileName.c_str(), NcFile::Replace);
589       // check if writing was successful
590       if (!dataFile.is_valid())
591            throw DataException("Error - DataTagged:: opening of netCDF file for output failed.");
592       if (!dataFile.add_att("type_id",1) )
593            throw DataException("Error - DataTagged:: appending data type to netCDF file failed.");
594       if (!dataFile.add_att("rank",rank) )
595            throw DataException("Error - DataTagged:: appending rank attribute to netCDF file failed.");
596       if (!dataFile.add_att("function_space_type",type))
597            throw DataException("Error - DataTagged:: appending function space attribute to netCDF file failed.");
598       ndims=rank+1;
599       if ( rank >0 ) {
600           dims[0]=shape[0];
601           if (! (ncdims[0] = dataFile.add_dim("d0",shape[0])) )
602                throw DataException("Error - DataTagged:: appending ncdimsion 0 to netCDF file failed.");
603       }
604       if ( rank >1 ) {
605           dims[1]=shape[1];
606           if (! (ncdims[1] = dataFile.add_dim("d1",shape[1])) )
607                throw DataException("Error - DataTagged:: appending ncdimsion 1 to netCDF file failed.");
608       }
609       if ( rank >2 ) {
610           dims[2]=shape[2];
611           if (! (ncdims[2] = dataFile.add_dim("d2", shape[2])) )
612                throw DataException("Error - DataTagged:: appending ncdimsion 2 to netCDF file failed.");
613       }
614       if ( rank >3 ) {
615           dims[3]=shape[3];
616           if (! (ncdims[3] = dataFile.add_dim("d3", shape[3])) )
617                throw DataException("Error - DataTagged:: appending ncdimsion 3 to netCDF file failed.");
618       }
619       const DataTagged::DataMapType& thisLookup=getTagLookup();
620       DataTagged::DataMapType::const_iterator i;
621       DataTagged::DataMapType::const_iterator thisLookupEnd=thisLookup.end();
622       int ntags=1;
623       for (i=thisLookup.begin();i!=thisLookupEnd;i++) ntags++;
624       int* tags =(int*) esysUtils::malloc(ntags*sizeof(int));
625       int c=1;
626       tags[0]=-1;
627       for (i=thisLookup.begin();i!=thisLookupEnd;i++) tags[c++]=i->first;
628       dims[rank]=ntags;
629       if (! (ncdims[rank] = dataFile.add_dim("num_tags", dims[rank])) )
630       {
631           esysUtils::free(tags);
632               throw DataException("Error - DataTagged:: appending num_tags to netCDF file failed.");
633       }
634       if (! ( tags_var = dataFile.add_var("tags", ncInt, ncdims[rank])) )
635       {
636        esysUtils::free(tags);
637            throw DataException("Error - DataTagged:: appending tags to netCDF file failed.");
638       }
639       if (! (tags_var->put(tags,dims[rank])) )
640       {
641        esysUtils::free(tags);
642            throw DataException("Error - DataTagged:: copy tags to netCDF buffer failed.");
643       }
644       if (! ( var = dataFile.add_var("data", ncDouble, ndims, ncdims)) )
645       {
646        esysUtils::free(tags);
647            throw DataException("Error - DataTagged:: appending variable to netCDF file failed.");
648       }
649       if (! (var->put(d_ptr,dims)) )
650       {
651        esysUtils::free(tags);
652            throw DataException("Error - DataTagged:: copy data to netCDF buffer failed.");
653       }
654       #else
655       throw DataException("Error - DataTagged:: dump is not configured with netCDF. Please contact your installation manager.");
656       #endif
657    }
658  }  // end of namespace  }  // end of namespace

Legend:
Removed from v.594  
changed lines
  Added in v.1628

  ViewVC Help
Powered by ViewVC 1.1.26