CoreFlows/examples/Python/MPI4PY/testTwoSimulations.py

   1 #!/usr/bin/env python3
   2 # -*-coding:utf-8 -*
   3
   4 #===============================================================================================================================
   5 # Name        : Tests of launching two independent simulations in parallel
   6 # Author      : Michaël Ndjinga
   7 # Copyright   : CEA Saclay 2021
   8 # Description :
   9 #================================================================================================================================
  10
  11 from mpi4py import MPI
  12 import numpy as np
  13 import solverlab
  14 from math import sin, pi
  15
  16 def StationaryDiffusionEquation_2DEF_StructuredTriangles_par(split_direction, rank):
  17         spaceDim = 2;
  18         # Prepare for the mesh
  19         print("Processor ", rank, " : Building mesh " );
  20         xinf = 0 ;
  21         xsup=1.0;
  22         yinf=0.0;
  23         ysup=1.0;
  24         nx=20;
  25         ny=20;
  26         M=solverlab.Mesh(xinf,xsup,nx,yinf,ysup,ny,split_direction)#Regular triangular mesh
  27         # set the limit field for each boundary
  28         eps=1e-6;
  29         M.setGroupAtPlan(xsup,0,eps,"Bord1")
  30         M.setGroupAtPlan(xinf,0,eps,"Bord2")
  31         M.setGroupAtPlan(ysup,1,eps,"Bord3")
  32         M.setGroupAtPlan(yinf,1,eps,"Bord4")
  33
  34         print("Processor ", rank, " : Built a regular triangular 2D mesh from a square mesh with ", nx,"x" ,ny, " cells.")
  35         print("Processor ", rank, " : Each square was split in two in direction ",split_direction)
  36         FEComputation=True
  37         myProblem = solverlab.StationaryDiffusionEquation(spaceDim,FEComputation);
  38         myProblem.setMesh(M);
  39
  40     # set the limit value for each boundary
  41         T1=0;
  42         T2=0;
  43         T3=0;
  44         T4=0;
  45
  46         myProblem.setDirichletBoundaryCondition("Bord1",T1)
  47         myProblem.setDirichletBoundaryCondition("Bord2",T2)
  48         myProblem.setDirichletBoundaryCondition("Bord3",T3)
  49         myProblem.setDirichletBoundaryCondition("Bord4",T4)
  50
  51         #Set the right hand side function
  52         my_RHSfield = solverlab.Field("RHS_field", solverlab.NODES, M, 1)
  53         for i in range(M.getNumberOfNodes()):
  54                 Ni= M.getNode(i)
  55                 x = Ni.x()
  56                 y = Ni.y()
  57
  58                 my_RHSfield[i]=2*pi*pi*sin(pi*x)*sin(pi*y)#mettre la fonction definie au second membre de l'edp
  59
  60         myProblem.setHeatPowerField(my_RHSfield)
  61         myProblem.setLinearSolver(solverlab.GMRES,solverlab.ILU);
  62
  63     # name of result file
  64         fileName = "StationnaryDiffusion_2DEF_StructuredTriangles"+str(rank);
  65
  66     # computation parameters
  67         myProblem.setFileName(fileName);
  68
  69     # Run the computation
  70         myProblem.initialize();
  71         print("Processor ", rank, " : Running python "+ fileName );
  72
  73         ok = myProblem.solveStationaryProblem();
  74         if (not ok):
  75                 print( "Python simulation of " + fileName + "  failed ! " );
  76                 pass
  77         else:
  78                 ####################### Postprocessing #########################
  79                 my_ResultField = myProblem.getOutputTemperatureField()
  80                 #The following formulas use the fact that the exact solution is equal the right hand side divided by 2*pi*pi
  81                 max_abs_sol_exacte=max(my_RHSfield.max(),-my_RHSfield.min())/(2*pi*pi)
  82                 max_sol_num=my_ResultField.max()
  83                 min_sol_num=my_ResultField.min()
  84                 erreur_abs=0
  85                 for i in range(M.getNumberOfNodes()) :
  86                         if erreur_abs < abs(my_RHSfield[i]/(2*pi*pi) - my_ResultField[i]) :
  87                                 erreur_abs = abs(my_RHSfield[i]/(2*pi*pi) - my_ResultField[i])
  88
  89                 print("Processor ", rank, " : Absolute error = max(| exact solution - numerical solution |) = ",erreur_abs )
  90                 print("Processor ", rank, " : Relative error = max(| exact solution - numerical solution |)/max(| exact solution |) = ",erreur_abs/max_abs_sol_exacte)
  91                 print("Processor ", rank, " : Maximum numerical solution = ", max_sol_num, " Minimum numerical solution = ", min_sol_num)
  92
  93                 assert erreur_abs/max_abs_sol_exacte <1.
  94                 pass
  95
  96         print("Processor ", rank, " : ------------ !!! End of calculation !!! -----------" );
  97
  98         myProblem.terminate();
  99         return erreur_abs/max_abs_sol_exacte
 100
 101 if __name__ == """__main__""":
 102         comm = MPI.COMM_WORLD
 103         size = comm.Get_size()
 104         rank = comm.Get_rank()
 105
 106         if(size!=2):
 107                 raise ValueError("Processor ", rank, " : aborting.\n Simulation should done on two processors.\n", size, " processors given")
 108
 109         print("My rank is ", rank, " among ", size, "processors ")
 110
 111         my_relative_error=StationaryDiffusionEquation_2DEF_StructuredTriangles_par(rank, rank)
 112
 113         if rank == 0:
 114             comm.send(my_relative_error, dest=1, tag=11)
 115             other_relative_error = comm.recv(source=1, tag=17)
 116         elif rank == 1:
 117             other_relative_error = comm.recv(source=0, tag=11)
 118             comm.send(my_relative_error, dest=0, tag=17)
 119
 120         print("Processor ", rank, " : Difference between the two processor relative errors is ", abs(my_relative_error-other_relative_error) )
 121         #print("Processor ", rank, " : Difference between the two processors is ", (my_ResultField-other_ResultField).normMax()[0] )
 122