doc/tutorial/medcoupling_dataarray1_en.rst

   1
   2 Manipulate DataArray objects
   3 ----------------------------
   4
   5 DataArrays (DataArrayInt and DataArrayDouble) are used in MEDCoupling to store contiguous arrays in memory. Each array entry is a tuple of data.
   6 They form the base of many data processing performed by MEDCoupling. It is thus useful to be able to manipulate them with ease.
   7
   8 DataArrayDouble are often used to manipulate fields in an optimized fashion, as will be shown later on.
   9
  10 Objective
  11 ~~~~~~~~~
  12
  13 In this exercise the aim is to create a set of 2D Cartesian coordinates. The first component will be called X (with unit "m") and the second Y (same unit). The coordinates represent 7 regular hexagons (all inscribed in a 3cm radius circle).
  14 The 7 hexagons are built from a centered reference hexagon (3.4; 4.4).
  15
  16 Around this reference, 6 translated copies are created and each copy shares exactly one edge with the reference.
  17 The common nodes (tuples) are then merged. This shows how to use indirections, concept widely used in unstructured meshes.
  18
  19 .. image:: images/DataArrayDouble_1.jpg
  20
  21 Concepts at hand:
  22
  23 * Create an instance of DataArrayDouble
  24 * Display a DataArrayDouble instance and invoke the getValue() method to convert it to a list
  25 * Use convenient notations da[...]
  26 * Learn about renumbering (old2new convention)
  27 * Invoke services from findCommonTuples
  28
  29 Starting the implementation
  30 ~~~~~~~~~~~~~~~~~~~~~~~~~~~
  31
  32 To start with, import the whole python module medcoupling and the standard math module. ::
  33
  34         from medcoupling import *
  35         import math
  36
  37 This makes the following available:
  38
  39 * all MEDCoupling classes
  40 * all enum (ON_CELLS, ON_NODES, ONE_TIME...)
  41 * all static methods
  42
  43 Create a DataArrayDouble instance containing 6 tuples
  44 ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
  45
  46 The target here is to create one DataArrayDouble containing the complete set of coordinates for one regular hexagon. ::
  47
  48         d=DataArrayDouble(6,2)
  49
  50 Which is equivalent to ::
  51
  52         d=DataArrayDouble()
  53         d.alloc(6,2)
  54
  55 or to ::
  56
  57         d=DataArrayDouble(12)
  58         d.rearrange(2)
  59
  60 .. note:: d contains 12 values grouped in 6 tuples containing each 2 components.
  61           The values in d are still un-initialized.
  62
  63 Initialize an instance of DataArrayDouble
  64 ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
  65
  66 Assign the value 3.0 (radius) to the first component of all tuples. ::
  67
  68         d[:,0]=3.
  69
  70 Initialize the 2nd component of each tuple i with the index value i. ::
  71
  72         d[:,1]=list(range(6))
  73
  74 Multiply the 2nd component by pi/3. ::
  75
  76         d[:,1]*=math.pi/3.
  77
  78 .. note:: d now contains the polar coordinates of our regular hexagon centered at 0,0.
  79
  80 Convert from polar to Cartesian by invoking fromPolarToCart(), and assign the result to d. ::
  81
  82         d=d.fromPolarToCart()
  83
  84 .. note:: fromPolarToCart() generates a new instance and d now points to it.
  85
  86 Assign the correct component information of d.::
  87
  88         d.setInfoOnComponents(["X [m]","Y [m]"])
  89
  90 Display the values only in form of a Python list. ::
  91
  92         print(d.getValues())
  93
  94 Display d. ::
  95
  96         print(d)
  97
  98 Verify that for each tuple of d its norm (returned by the magnitude() method) is indeed 3 (with an error less than 1.e-12):
  99 ::
 100
 101         print(d.magnitude().isUniform(3.,1e-12))
 102
 103
 104 Duplication and aggregation of DataArrayDouble
 105 ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
 106
 107 translationToPerform contains a list of vectors of size 2. This list has length 7 (7 hexagons) Cette liste de taille 7 (7 hexagones) and describes the appropriate translation for each copy of d.
 108
 109 Execute directly the following lines ::
 110
 111         radius=3.
 112         translationToPerform=[[0.,0.],[3./2.*radius,-radius*math.sqrt(3.)/2],[3./2.*radius,radius*math.sqrt(3.)/2],[0.,radius*math.sqrt(3.)],[-3./2.*radius,radius*math.sqrt(3.)/2],[-3./2.*radius,-radius*math.sqrt(3.)/2],[0.,-radius*math.sqrt(3.)]]
 113
 114
 115 Create the len(translationToPerform) copies of d and apply the corresponding translation.  ::
 116
 117         ds=len(translationToPerform)*[None]
 118         for pos,t in enumerate(translationToPerform):
 119           ds[pos]=d[:]
 120           ds[pos]+=t
 121           pass
 122
 123 An alternative (and more compact) way to do it : ::
 124
 125         ds=[d.deepCopy() for i in list(range(len(translationToPerform)))]
 126         for (elt,t) in zip(ds,translationToPerform) : elt+=t
 127
 128 Aggregating DataArrayDouble
 129 ~~~~~~~~~~~~~~~~~~~~~~~~~~~
 130 From the list of DataArrayDouble instances 'ds', construct the aggregation 'd2' (a new instance of DataArrayDouble).
 131 ::
 132
 133         d2=DataArrayDouble.Aggregate(ds)
 134
 135 .. note:: d2 now contains the full tuple set (6*7=42 tuples of 2 components each) which was in ds with the same order. This might sound obvious but the aggregation of meshes and fields is based on the exact same principle thus facilitating data access and indexing. This is a key difference to the MED file model, as will be seen later on.
 136
 137 Find identical tuples in d2
 138 ~~~~~~~~~~~~~~~~~~~~~~~~~~~
 139
 140 d2 contains 42 tuples but some of them are duplicated. To identify those (with a 1e-12 precision) invoke DataArrayDouble.findCommonTuples().
 141 Use help(DataArrayDouble.findCommonTuples) to show the prototype of the method and store the result in c and cI::
 142
 143         oldNbOfTuples=d2.getNumberOfTuples()
 144         c,cI=d2.findCommonTuples(1e-12)
 145
 146 c is an array containing all the common nodes grouped sequentially. cI contains the pointers to c allowing to identify all the nodes in a group 'i' (with 'i' in [0, 12) ). Thus the tuple ids of group 'i' start at index cI[i] in the list c, and finish at index cI[i+1].
 147
 148 .. note:: DataArrayDouble.findCommonTuples() thus returns 2 values: an array containing a list of common tuples, and an index array which is used to navigate in the first list. This is a very classical return signature in MEDCoupling and often appears when manipulating unstructured meshes.
 149
 150 Manipulate DataArrayInt couple representing a pair (Data,IndexData)
 151 ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
 152
 153 The number of common tuples is len(cI)-1, i.e. 12 in our case.
 154 Get the list of tuple ids forming the first group of common nodes and store the result in tmp. ::
 155
 156         tmp=c[cI[0]:cI[0+1]]
 157         print(tmp)
 158
 159 Check the result: all the tuples stored in tmp point to identical coordinates in d2.::
 160
 161         print(d2[tmp])
 162
 163 .. note:: we see the tuple (3.,0.) repeated 3 times (with an error margin below 1e-12).
 164
 165 Now we will deduce from the 3 variables oldNbOfTuples, c and cI the number of truly different tuples in d2.
 166 To this end we compute the number of repetitions in d2 and we subtract this from oldNbOfTuples.
 167
 168 To get the repetition count invoke the method DataArrayInt.deltaShiftIndex(). It returns the size of each group.
 169 Store the result in a. ::
 170
 171         a=cI.deltaShiftIndex()
 172
 173 The number of repetitions in d2 for each group can now be stored in b. ::
 174
 175         b=a-1
 176
 177 Finally the new count of tuples can be computed as follows. ::
 178
 179         myNewNbOfTuples=oldNbOfTuples-sum(b.getValues())
 180
 181 Build old to new array from c and cI
 182 ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
 183
 184 c and cI define a surjective function from a starting set X containing 42 tuples (oldNbOfTuples) to a set Y containing 24 tuples (myNewNbOfTuples).
 185
 186 .. image:: images/SurjectionDataArray.png
 187
 188 This surjective function may also be represented thanks to an "old-2-new" array.
 189 This storage method takes the form of a DataArrayInt 'o2n' made of Card(X) tuples (in our case 42) with one component each. For each tuple (element) indexed by i in o2n, o2n[i] contains the new tuple id in Y.
 190
 191 The format 'old-2-new' is systematically used for all renumbering operations (one-to-one correspondence).
 192
 193 The static method DataArrayInt.ConvertIndexArrayToO2N() performs the conversion from one storage mode to the other (c, cI to o2n).
 194 We get for free the number of elements in Y, i.e. the variable newNbOfTuples. ::
 195
 196         o2n,newNbOfTuples=DataArrayInt.ConvertIndexArrayToO2N(oldNbOfTuples,c,cI)
 197         print("Have I got the right result? %s"%(str(myNewNbOfTuples==newNbOfTuples)))
 198
 199 Using o2n and newNbOfTuples invoke DataArrayDouble.renumberAndReduce() on d2. ::
 200
 201         d3=d2.renumberAndReduce(o2n,newNbOfTuples)
 202
 203 This method has one drawback: we don't know for each group of common tuples in d2 what id was finally retained.
 204 For example: in the group 0 we know that tuples 0, 8 and 16 are identical (tmp.getValues()), but we don't know which of the three has been retained in d3.
 205
 206 To make this choice explicit we use the new-2-old format. This storage mode is represented by a DataArrayInt n2o containing Card(Y) tuples (in our case 24) with one component each.
 207 For each tuple (element) with id 'i' in n2o, n2o[i] contains the tuple id which was chosen in X.
 208
 209 The method DataArrayInt.invertArrayO2N2N2O() allows to switch between the two representations.
 210 Try it on the variable o2n. ::
 211
 212         n2o=o2n.invertArrayO2N2N2O(newNbOfTuples)
 213
 214 Using n2o we can deduce d3_bis from d2. ::
 215
 216         d3_bis=d2[n2o]
 217         print("Have I got the right result (2)? %s"%(str(d3.isEqual(d3_bis,1e-12))))
 218
 219 Translate all tuples at once
 220 ~~~~~~~~~~~~~~~~~~~~~~~~~~~~
 221
 222 All tuples (or nodes) are to be translated by the vector [3.3,4.4]. ::
 223
 224         d3+=[3.3,4.4]
 225
 226 Build an unstructured mesh using d3 (coordinates) and o2n
 227 ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
 228
 229 Create an unstructured mesh m with coordinates d3. m has a mesh dimension of 2 ::
 230
 231         m=MEDCouplingUMesh("My7hexagons",2)
 232         m.setCoords(d3)
 233
 234 Now allocate the number of cells with an upper bound of the actual number of cells. ::
 235
 236         m.allocateCells(7)
 237
 238 Finally thanks to o2n we know the connectivity of all 7 hexagons using the coordinates stored in d3. ::
 239
 240         for i in list(range(7)):
 241           m.insertNextCell(NORM_POLYGON,o2n[6*i:6*(i+1)].getValues())
 242           pass
 243
 244 Check that m is coherent. ::
 245
 246          m.checkConsistencyLight()
 247
 248 To visually check m, write it in a VTU file ("My7hexagons.vtu") and display it in ParaVis. ::
 249
 250         m.writeVTK("My7hexagons.vtu")
 251
 252 .. note:: We write the data in VTU file, not in a MED file, because MEDCoupling doesn't use the MED file pre-requisite. Also note that although we output VTK data, MEDCoupling in itself doesn't depend on the VTK libraries.
 253
 254 Solution
 255 ~~~~~~~~
 256
 257 :ref:`python_testMEDCouplingdataarray1_solution`