Salome HOME
Indices are stored as mcIdType type instead of int to support switch to 64bits indexing
[tools/medcoupling.git] / src / ParaMEDMEM_Swig / test_InterpKernelDEC.py
index 27e4551de3d359d3e517d44a109b71181f571670..a65c4638cd0073da05336a9647f66e62cd8b446c 100755 (executable)
@@ -1,6 +1,6 @@
 #!/usr/bin/env python
 #  -*- coding: iso-8859-1 -*-
-# Copyright (C) 2007-2015  CEA/DEN, EDF R&D
+# Copyright (C) 2007-2019  CEA/DEN, EDF R&D
 #
 # This library is free software; you can redistribute it and/or
 # modify it under the terms of the GNU Lesser General Public
 #
 
 from ParaMEDMEM import *
+from MEDLoader import ReadUMeshFromFile
 import sys, os
 import unittest
 import math
+from mpi4py import MPI
+
 
 class ParaMEDMEMBasicsTest(unittest.TestCase):
     def testInterpKernelDEC_2D(self):
-        MPI_Init(sys.argv)
-        size = MPI_Comm_size(MPI_COMM_WORLD)
-        rank = MPI_Comm_rank(MPI_COMM_WORLD)
+        size = MPI.COMM_WORLD.size
+        rank = MPI.COMM_WORLD.rank
         if size != 5:
-            raise RuntimeError, "Expect MPI_COMM_WORLD size == 5"
-        print rank
+            raise RuntimeError("Expect MPI_COMM_WORLD size == 5")
+        print(rank)
         nproc_source = 3
-        procs_source = range( nproc_source )
-        procs_target = range( size - nproc_source + 1, size)
+        procs_source = list(range(nproc_source))
+        procs_target = list(range(size - nproc_source + 1, size))
 
         interface = CommInterface()
         target_group = MPIProcessorGroup(interface, procs_target)
@@ -56,35 +58,33 @@ class ParaMEDMEMBasicsTest(unittest.TestCase):
         filename_xml1 = os.path.join(data_dir, "share/resources/med/square1_split")
         filename_xml2 = os.path.join(data_dir, "share/resources/med/square2_split")
 
-        MPI_Barrier(MPI_COMM_WORLD)
+        MPI.COMM_WORLD.Barrier()
         if source_group.containsMyRank():
             filename = filename_xml1 + str(rank+1) + ".med"
             meshname = "Mesh_2_" + str(rank+1)
-            mesh=MEDLoader.ReadUMeshFromFile(filename,meshname,0)
+            mesh=ReadUMeshFromFile(filename,meshname,0)
             paramesh=ParaMESH(mesh,source_group,"source mesh")
             comptopo = ComponentTopology()
             parafield = ParaFIELD(ON_CELLS,NO_TIME,paramesh, comptopo)
-            parafield.getField().setNature(ConservativeVolumic)
+            parafield.getField().setNature(IntensiveMaximum)
             nb_local=mesh.getNumberOfCells()
             value = [1.0]*nb_local
             parafield.getField().setValues(value)
-            icocofield = ICoCoMEDField(mesh,parafield.getField())
-            dec.setMethod("P0")
+            icocofield = ICoCoMEDField(parafield.getField())
             dec.attachLocalField(icocofield)
             pass
         else:
             filename = filename_xml2 + str(rank - nproc_source + 1) + ".med"
             meshname = "Mesh_3_" + str(rank - nproc_source + 1)
-            mesh=MEDLoader.ReadUMeshFromFile(filename,meshname,0)
+            mesh=ReadUMeshFromFile(filename,meshname,0)
             paramesh=ParaMESH(mesh,target_group,"target mesh")
             comptopo = ComponentTopology()
             parafield = ParaFIELD(ON_CELLS,NO_TIME,paramesh, comptopo)
-            parafield.getField().setNature(ConservativeVolumic)
+            parafield.getField().setNature(IntensiveMaximum)
             nb_local=mesh.getNumberOfCells()
             value = [0.0]*nb_local
             parafield.getField().setValues(value)
-            icocofield = ICoCoMEDField(mesh,parafield.getField())
-            dec.setMethod("P0")
+            icocofield = ICoCoMEDField(parafield.getField())
             dec.attachLocalField(icocofield)
             pass
 
@@ -95,7 +95,7 @@ class ParaMEDMEMBasicsTest(unittest.TestCase):
             dec.sendData()
             dec.recvData()
             field_after_int=parafield.getVolumeIntegral(0,True);
-            self.failUnless(math.fabs(field_after_int-field_before_int)<1e-8)
+            self.assertTrue(math.fabs(field_after_int-field_before_int)<1e-8)
             pass
         else:
             dec.synchronize()
@@ -113,8 +113,8 @@ class ParaMEDMEMBasicsTest(unittest.TestCase):
         paramesh   =0
         parafield  =0
         icocofield =0
-        MPI_Barrier(MPI_COMM_WORLD)
-        MPI_Finalize()
+        MPI.COMM_WORLD.Barrier()
+        MPI.Finalize()
         pass
     pass