doc/tutorial/atestMEDCouplingRemapper.rst

   1
   2 .. _python_testMEDCouplingremapper1_solution:
   3
   4 Interpoler avec MEDCouplingRemapper
   5 ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
   6
   7 ::
   8
   9         import MEDCoupling as mc
  10         from MEDCouplingRemapper import MEDCouplingRemapper
  11         # Target mesh
  12         arr = mc.DataArrayDouble(11)
  13         arr.iota(0)
  14         trgMesh = mc.MEDCouplingCMesh()
  15         trgMesh.setCoords(arr,arr)
  16         trgMesh = trgMesh.buildUnstructured()
  17         # Source mesh
  18         arr = mc.DataArrayDouble(21)
  19         arr.iota(0)
  20         arr *= 0.5
  21         srcMesh = mc.MEDCouplingCMesh()
  22         srcMesh.setCoords(arr,arr)
  23         srcMesh = srcMesh.buildUnstructured()
  24         # Triangularize some cells in source
  25         tmp = srcMesh[:20]    # Extract a sub-part of srcMesh
  26         tmp.simplexize(0)
  27         srcMesh = mc.MEDCouplingUMesh.MergeUMeshes([tmp,srcMesh[20:]])
  28         # Interpolate cells to cells
  29         remap = MEDCouplingRemapper()
  30         remap.prepare(srcMesh,trgMesh,"P0P0")
  31         # Check matrix
  32         myMatrix = remap.getCrudeMatrix()
  33         print myMatrix
  34         sumByRows = mc.DataArrayDouble(len(myMatrix))
  35         for i,wIt in enumerate(sumByRows):
  36           su = 0.
  37           for it in myMatrix[i]:
  38             su += myMatrix[i][it]
  39           wIt[0] = su
  40         print "Is interpolation well prepared?", sumByRows.isUniform(1.,1e-12)
  41         # Source field construction
  42         srcField = mc.MEDCouplingFieldDouble(mc.ON_CELLS, mc.ONE_TIME)
  43         srcField.setMesh(srcMesh)
  44         srcField.fillFromAnalytic(1,"7-sqrt((x-5.)*(x-5.)+(y-5.)*(y-5.))")
  45         srcField.getArray().setInfoOnComponent(0, "powercell [W]")
  46         # Transfer field
  47         #remap.transferField(srcField, 1e300)
  48         srcField.setNature(mc.ConservativeVolumic)
  49         trgFieldCV = remap.transferField(srcField,1e300)
  50         # ConservativeVolumic
  51         integSource = srcField.integral(True)[0]
  52         integTarget =  trgFieldCV.integral(True)[0]
  53         print "ConservativeVolumic -- integrals: %lf == %lf" % (integSource, integTarget)
  54
  55         accSource = srcField.getArray().accumulate()[0]
  56         accTarget = trgFieldCV.getArray().accumulate()[0]
  57         print "ConservativeVolumic -- sums: %lf != %lf" % (accSource, accTarget)
  58         # IntegralGlobConstraint
  59         srcField.setNature(mc.IntegralGlobConstraint)
  60         trgFieldI = remap.transferField(srcField,1e300)
  61         #
  62         integSource = srcField.integral(True)[0]
  63         integTarget =  trgFieldI.integral(True)[0]
  64         print "IntegralGlobConstraint -- integrals: %lf != %lf" % (integSource, integTarget)
  65
  66         accSource = srcField.getArray().accumulate()[0]
  67         accTarget = trgFieldI.getArray().accumulate()[0]
  68         print "IntegralGlobConstraint -- sums: %lf == %lf" % (accSource, accTarget)