Back to index

salome-med  6.5.0
test_StructuredCoincidentDEC.py
Go to the documentation of this file.
00001 #!/usr/bin/env python
00002 #  -*- coding: iso-8859-1 -*-
00003 # Copyright (C) 2007-2012  CEA/DEN, EDF R&D
00004 #
00005 # This library is free software; you can redistribute it and/or
00006 # modify it under the terms of the GNU Lesser General Public
00007 # License as published by the Free Software Foundation; either
00008 # version 2.1 of the License.
00009 #
00010 # This library is distributed in the hope that it will be useful,
00011 # but WITHOUT ANY WARRANTY; without even the implied warranty of
00012 # MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
00013 # Lesser General Public License for more details.
00014 #
00015 # You should have received a copy of the GNU Lesser General Public
00016 # License along with this library; if not, write to the Free Software
00017 # Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307 USA
00018 #
00019 # See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
00020 #
00021 
00022 from ParaMEDMEM import *
00023 import sys, os
00024 import unittest
00025 import math
00026 
00027 class ParaMEDMEMBasicsTest2(unittest.TestCase):
00028     def testStructuredCoincidentDEC(self):
00029         MPI_Init(sys.argv)
00030         #
00031         size = MPI_Comm_size(MPI_COMM_WORLD)
00032         rank = MPI_Comm_rank(MPI_COMM_WORLD)
00033         #
00034         if size < 4:
00035             raise RuntimeError, "Expect MPI_COMM_WORLD size >= 4"
00036         #
00037         interface = CommInterface()
00038         #
00039         self_group   = MPIProcessorGroup(interface, rank, rank)
00040         target_group = MPIProcessorGroup(interface, 3, size-1)
00041         source_group = MPIProcessorGroup(interface, 0, 2)
00042         #
00043         mesh      = 0
00044         support   = 0
00045         paramesh  = 0
00046         parafield = 0
00047         comptopo  = 0
00048         icocofield= 0
00049         #
00050         data_dir = os.environ['MED_ROOT_DIR']
00051         tmp_dir  = os.environ['TMP']
00052         if tmp_dir == '':
00053             tmp_dir = "/tmp"
00054             pass
00055         
00056         filename_xml1    = data_dir + "/share/salome/resources/med/square1_split"
00057         filename_2       = data_dir + "/share/salome/resources/med/square1.med"
00058         filename_seq_wr  = tmp_dir + "/"
00059         filename_seq_med = tmp_dir + "/myWrField_seq_pointe221.med"
00060         
00061         dec = StructuredCoincidentDEC(source_group, target_group)
00062         MPI_Barrier(MPI_COMM_WORLD)
00063         if source_group.containsMyRank():
00064             filename = filename_xml1 + str(rank+1) + ".med"
00065             meshname = "Mesh_2_" + str(rank+1)
00066             mesh=MEDLoader.ReadUMeshFromFile(filename,meshname,0)
00067             paramesh=ParaMESH(mesh,source_group,"source mesh")
00068             comptopo=ComponentTopology(6)
00069             parafield=ParaFIELD(ON_CELLS,NO_TIME,paramesh,comptopo)
00070             parafield.getField().setNature(ConservativeVolumic)
00071             nb_local=mesh.getNumberOfCells()
00072             global_numbering=paramesh.getGlobalNumberingCell2()
00073             value = []
00074             for ielem in range(nb_local):
00075                 for icomp in range(6):
00076                     value.append(global_numbering[ielem]*6.0+icomp);
00077                     pass
00078                 pass
00079             parafield.getField().setValues(value)
00080             icocofield = ICoCoMEDField(mesh,parafield.getField())
00081             dec.setMethod("P0")  
00082             dec.attachLocalField(parafield)      
00083             dec.synchronize()
00084             dec.sendData()
00085             pass
00086 
00087         if target_group.containsMyRank():
00088             meshname2 = "Mesh_2"
00089             mesh=MEDLoader.ReadUMeshFromFile(filename_2, meshname2,0)
00090             paramesh=ParaMESH(mesh, self_group, "target mesh")
00091             comptopo=ComponentTopology(6,target_group)
00092             parafield=ParaFIELD(ON_CELLS,NO_TIME,paramesh, comptopo)
00093             parafield.getField().setNature(ConservativeVolumic)
00094             nb_local=mesh.getNumberOfCells()
00095             value = [0.0]*(nb_local*comptopo.nbLocalComponents())
00096             parafield.getField().setValues(value)
00097             icocofield = ICoCoMEDField(mesh,parafield.getField())
00098             dec.setMethod("P0")
00099             dec.attachLocalField(parafield)
00100             dec.synchronize()
00101             dec.recvData()
00102             recv_value = parafield.getField().getArray().getValues()
00103             for i in range(nb_local):
00104                 first=comptopo.firstLocalComponent()
00105                 for icomp in range(comptopo.nbLocalComponents()):
00106                     self.failUnless(math.fabs(recv_value[i*comptopo.nbLocalComponents()+icomp]-
00107                                               (float)(i*6+icomp+first))<1e-12)
00108                     pass
00109                 pass
00110             pass
00111         comptopo=0
00112         interface = 0
00113         mesh       =0
00114         support    =0
00115         paramesh   =0
00116         parafield  =0
00117         icocofield =0
00118         dec=0
00119         self_group =0
00120         target_group = 0
00121         source_group = 0
00122         MPI_Barrier(MPI_COMM_WORLD)
00123         MPI_Finalize()
00124         print "End of test StructuredCoincidentDEC"
00125         pass
00126 
00127     
00128 unittest.main()