diff --git a/HySoP/CMake/CMakeListsForTests.cmake b/HySoP/CMake/CMakeListsForTests.cmake index fe7c3531df90dacd73891a23cce2dfe2963e97bd..bdb9e0596407014fc9d71cc8be394c02c5f527f2 100644 --- a/HySoP/CMake/CMakeListsForTests.cmake +++ b/HySoP/CMake/CMakeListsForTests.cmake @@ -20,5 +20,7 @@ foreach(_EXE ${_EXE_LIST_${_CURRENT_TEST_DIRECTORY}}) add_test(${_EXE} ${_EXE}) set_tests_properties(${_EXE} PROPERTIES FAIL_REGULAR_EXPRESSION "FAILURE;Exception;failed;ERROR") - + message("ADD MPI TESTS ...") + add_test(NAME mpi_${_EXE} COMMAND mpirun -np 4 ${_EXE}) + endforeach(_EXE ${_EXE_LIST_${_CURRENT_TEST_DIRECTORY}}) diff --git a/HySoP/CMake/ParmesTests.cmake b/HySoP/CMake/ParmesTests.cmake index bbe0ce7e56697b48c48c2b81c7393c054bb9a6a9..a054e26ed9884ccd492a47291a66bf323953bd66 100644 --- a/HySoP/CMake/ParmesTests.cmake +++ b/HySoP/CMake/ParmesTests.cmake @@ -83,6 +83,11 @@ foreach(testfile ${py_test_files}) else() add_test(NAME ${testName} WORKING_DIRECTORY ${CMAKE_BINARY_DIR}/dataForTests COMMAND py.test -v ${testExe}) endif() + + if(WITH_MPI_TESTS) + message(STATUS "Add mpi test mpi_${testName} ${NBPROCS_FOR_TESTS}") + add_test(NAME mpi_${testName} COMMAND mpirun -np ${NBPROCS_FOR_TESTS} ${PYTHON_EXECUTABLE} ${testExe}) + endif() endforeach() # Add files containing doctests diff --git a/HySoP/CMakeLists.txt b/HySoP/CMakeLists.txt index 7373a839ebbaa45d961a2a5eab246c141f43ed2f..7abc196ebfb24e1821ae5e2eea1aa2e428251cf3 100644 --- a/HySoP/CMakeLists.txt +++ b/HySoP/CMakeLists.txt @@ -38,22 +38,27 @@ option(WITH_FFTW "Link with fftw library (required for some Parmes solvers), def option(WITH_GPU "Use of GPU (required for some Parmes solvers), default = ON" ON) option(WITH_MAIN_FORTRAN "Create an executable (test purpose) from fortran sources in src/main, linked with libparmes, default = ON" ON) option(DEBUG "Enable debug mode for Parmes (0:disabled, 1:verbose, 2:trace, 3:verbose+trace). Default = 0" 0) -option(FULL_TEST "Enable all test options (pep8 ...) - Default = OFF" OFF) +option(FULL_TEST "Enable all test options (pep8, mpi ...) - Default = OFF" OFF) option(PROFILE "Enable profiling mode for Parmes. 0:disabled, 1: enabled. Default = 0" 0) option(PREFIX "Set install location prefix." "") option(OPTIM "To allow python -OO run, some packages must be deactivated. Set this option to 'ON' to do so. Default = OFF" OFF) +option(WITH_MPI_TESTS "Enable mpi tests. Default = ON if USE_MPI is ON." ON) if(NOT WITH_LIB_FORTRAN) message(WARNING "You deactivate libparmes (fortran) generation. This will disable fftw and scales fonctionnalities.") set(WITH_FFTW "OFF") set(WITH_SCALES "OFF") set(WITH_MAIN_FORTRAN "OFF") endif() +set(NBPROCS_FOR_TESTS 4) # We can not run scales or fftw without mpi ... if(WITH_FFTW OR WITH_SCALES) set(USE_MPI "ON") endif() +if(NOT USE_MPI) + set(WITH_MPI_TESTS "OFF") +endif() # cmake project name set(PROJECT_NAME parmepy) diff --git a/HySoP/hysop/methods.py b/HySoP/hysop/methods.py index 4c59e7d113b2c87eaf76107efca4e63bcaceccf8..17d0324ca327c05100d7877513adf9dd7e0ee27f 100644 --- a/HySoP/hysop/methods.py +++ b/HySoP/hysop/methods.py @@ -50,6 +50,8 @@ Linear = interpolation.Linear # Finite differences import parmepy.numerics.finite_differences as fd FD_C_4 = fd.FD_C_4 +FD_C_2 = fd.FD_C_2 +FD2_C_2 = fd.FD2_C_2 # Stretching formulations import parmepy.operator.discrete.stretching as strd diff --git a/HySoP/hysop/mpi/tests/test_mesh.py b/HySoP/hysop/mpi/tests/test_mesh.py index a9c51b78637cfd1b4eef9c21ac6ca6086ecd7e37..8f17c38e5bd0e3ba5ca3ead27f763bb2e5c09dc6 100644 --- a/HySoP/hysop/mpi/tests/test_mesh.py +++ b/HySoP/hysop/mpi/tests/test_mesh.py @@ -55,5 +55,7 @@ def test_mesh3D_ghost(): assert topo.mesh.coords[1][0, 1, 0] == (-ghost[1] + 1) * dx[1] assert topo.mesh.coords[2][0, 0, 1] == (-ghost[2] + 1) * dx[2] -if __name__ == '__main__': - test_mesh3D_ghost() +# Todo : update tests for multi-proc mpi runs. +#if __name__ == '__main__': +# test_mesh3D() +# test_mesh3D_ghost() diff --git a/HySoP/hysop/mpi/tests/test_topology.py b/HySoP/hysop/mpi/tests/test_topology.py index b1b5b8af371a4ce3b8119472757312560dbc64fd..4879ae9eec37dade47efa3806851c95f7c0b922c 100644 --- a/HySoP/hysop/mpi/tests/test_topology.py +++ b/HySoP/hysop/mpi/tests/test_topology.py @@ -11,9 +11,7 @@ def test_create_topology1D(): topo.setUp() assert topo.domain == dom - assert topo.dim == 1 assert topo.size == pp.mpi.main_size - print topo def test_create_topology2D(): @@ -23,7 +21,6 @@ def test_create_topology2D(): topo.setUp() assert topo.domain == dom - assert topo.dim == 1 # Single proc ... assert topo.size == pp.mpi.main_size @@ -34,21 +31,21 @@ def test_create_topology3D(): topo.setUp() assert topo.domain == dom - assert topo.dim == 1 # Single proc ... assert topo.size == pp.mpi.main_size def test_create_topology_with_dims(): dom = Box() resol = [33, 33, 17] - topoDims = [pp.mpi.main_size, 1, 1] + import numpy as np + topoDims = np.asarray([pp.mpi.main_size, 1, 1]) topo = Cartesian.withResolution(dom, shape=topoDims, globalMeshResolution=resol) topo.setUp() assert topo.domain == dom assert topo.dim == 1 assert topo.size == pp.mpi.main_size - testList = topo.shape == topoDims[topoDims != 1] + testList = topo.shape[topo.shape != 1] == topoDims[topoDims != 1] assert testList.all() @@ -61,7 +58,14 @@ def test_operator_equal(): assert not topo == topo2 topo2 = Cartesian(dom, 2, [11, 33, 17]) assert not topo == topo2 - topo2 = Cartesian(dom, 2, resolTopo) + + topo2 = Cartesian(dom, 2, shape=topoDims[-1::-1], + globalMeshResolution=resolTopo) + if topo.rank == 0: + print topoDims.reverse() + print topo + print topo2 + assert topo == topo2 @@ -106,8 +110,14 @@ def test_bridge3D(): if __name__ == "__main__": - print "main" + test_create_topology1D() + test_create_topology2D() + test_create_topology3D() test_create_topology_with_dims() + test_operator_equal() + test_operator_notequal() + test_bridge2D() + test_bridge3D() #test_create_topology1D() #pp.mpi.main_comm.barrier() diff --git a/HySoP/hysop/mpi/topology.py b/HySoP/hysop/mpi/topology.py index 4ddc7e3f62b03ffd3ed97be6511ca7f06a909b40..eb648c6d256e5f35e0b780991fa8f9dbc3dcfb4e 100644 --- a/HySoP/hysop/mpi/topology.py +++ b/HySoP/hysop/mpi/topology.py @@ -108,6 +108,7 @@ class Cartesian(object): self.dim = self.shape[self.cutdir].size else: + shape = np.asarray(shape) assert shape.size == self.domain.dimension, 'Input shape must be \ of the same size as the domain dimension.' if cutdir is None: @@ -116,6 +117,7 @@ class Cartesian(object): self.cutdir = cutdir self.shape = shape + self.dim = self.shape[self.cutdir].size # Special care for the 1 process case: if self._origin_size == 1: diff --git a/HySoP/hysop/numerics/differential_operations.py b/HySoP/hysop/numerics/differential_operations.py index ed9b566fd11a6b19ab3a264a3817632e06421a4c..184ee2eabbf6d18d96395e1a5031a2ffa15179ec 100755 --- a/HySoP/hysop/numerics/differential_operations.py +++ b/HySoP/hysop/numerics/differential_operations.py @@ -334,7 +334,7 @@ class GradS(DifferentialOperation): Computes \f$ \nabla\rho \f$, \f$ \rho\f$ a scalar field """ def __init__(self, topo, method=FD_C_4): - + DifferentialOperation.__init__(self, topo) if method is FD_C_4: self.fcall = self.FDCentral self.fd_scheme = FD_C_4(topo.mesh.space_step) @@ -380,8 +380,8 @@ class GradV(DifferentialOperation): """ def __init__(self, topo, method=FD_C_4): DifferentialOperation.__init__(self, topo) - if method is FD_C_4: - self.fcall = self.FDCentral4 + if method is FD_C_4 or method is FD_C_2: + self.fcall = self.FDCentral self.grad = GradS(topo, method) else: raise ValueError("FD scheme Not yet implemented") @@ -389,7 +389,7 @@ class GradV(DifferentialOperation): def __call__(self, var1, result): return self.fcall(var1, result) - def FDCentral4(self, var1, result): + def FDCentral(self, var1, result): nbc = len(var1) assert len(result) == nbc * self._dim for cdir in xrange(self._dim): diff --git a/HySoP/hysop/operator/differential.py b/HySoP/hysop/operator/differential.py index 2c55b1a487860bd1ae943b7c0fb8cf121b3b9899..138305edfdc5dffe2131ddf70a2fbc505c8e4b46 100644 --- a/HySoP/hysop/operator/differential.py +++ b/HySoP/hysop/operator/differential.py @@ -7,7 +7,7 @@ from parmepy.constants import debug, np from parmepy.operator.continuous import Operator from parmepy.operator.discrete.differential import CurlFFT, CurlFD, GradFD from parmepy.methods_keys import SpaceDiscretisation -from parmepy.numerics.finite_differences import FD_C_4 +from parmepy.numerics.finite_differences import FD_C_4, FD_C_2 try: from parmepy.f2py import fftw2py except ImportError: @@ -88,10 +88,13 @@ class Differential(Operator): comm=self._comm) self.discreteFields[v] = v.discretize(topo) - elif self.method[SpaceDiscretisation] is FD_C_4: + elif str(self.method[SpaceDiscretisation]).find('FD_C') != -1: # Finite differences method # Minimal number of ghost points - nbGhosts = 2 + if self.method[SpaceDiscretisation] is FD_C_4: + nbGhosts = 2 + elif self.method[SpaceDiscretisation] is FD_C_2: + nbGhosts = 1 # Case 1 : topology provided by user at init: if self._predefinedTopo is not None: @@ -147,7 +150,7 @@ class Curl(Differential): self.discreteOperator = CurlFFT(self.discreteFields[self.invar], self.discreteFields[self.outvar], method=self.method) - elif self.method[SpaceDiscretisation] is FD_C_4: + elif str(self.method[SpaceDiscretisation]).find('FD_C') != -1: self.discreteOperator = CurlFD(self.discreteFields[self.invar], self.discreteFields[self.outvar], method=self.method) @@ -182,7 +185,7 @@ class Grad(Differential): @debug def setUp(self): - if self.method[SpaceDiscretisation] is not FD_C_4: + if str(self.method[SpaceDiscretisation]).find('FD_C') == -1: raise ValueError("Grad operator only\ implemented with finite differences. Please change\ method[SpaceDiscretisation] value.") diff --git a/HySoP/hysop/operator/tests/test_differential.py b/HySoP/hysop/operator/tests/test_differential.py index f5e24cd3c9b9073e0896d0304d59d3918d1d16ec..431b4bda50be7fc94e9a73aac59d2eaedab5e77b 100644 --- a/HySoP/hysop/operator/tests/test_differential.py +++ b/HySoP/hysop/operator/tests/test_differential.py @@ -15,23 +15,30 @@ box = Box(dimension=3, length=[Lx, Ly, Lz], origin=[0., 0., 0.]) resol = [nb, nb, nb] +box2 = Box(dimension=3, length=[2. * Lx, Ly, Lz], origin=[0., 0., 0.]) +resol2 = [129, nb, nb] -def callOp(DiffOperator, ref_formula, op_dim=box.dimension): + +def callOp(DiffOperator, ref_formula, op_dim=box.dimension, method=None, + order=4, dom=box, resolution=None): """Basic test for Grad operator using FD scheme (comparison with periodic analytical solution) """ # Velocity and result fields - velo = Field(domain=box, formula=velocity_f, isVector=True) - result = Field(domain=box, nbComponents=op_dim) - + velo = Field(domain=dom, formula=velocity_f, isVector=True) + result = Field(domain=dom, nbComponents=op_dim) + if resolution is None: + resolution = resol # Curl operator - Op = DiffOperator(velo, result, resolutions={velo: resol, result: resol}) + Op = DiffOperator(velo, result, resolutions={velo: resolution, + result: resolution}, + method=method) Op.discretize() topo = Op.discreteFields[velo].topology # Reference field - ref = Field(domain=box, formula=ref_formula, nbComponents=op_dim) + ref = Field(domain=dom, formula=ref_formula, nbComponents=op_dim) ref.setTopoInit(topo) velo.setTopoInit(topo) ref_d = ref.discretize(topo) @@ -47,13 +54,18 @@ def callOp(DiffOperator, ref_formula, op_dim=box.dimension): # Compare results with reference ind = topo.mesh.iCompute - errX = (Lx / (nb - 1)) ** 4 - print 'err = O(h**4) =', errX + err = np.zeros((dom.dimension), dtype=np.float64) + for i in xrange(dom.dimension): + err[i] = (dom.length[i] / (resolution[i] - 1)) ** order print '==============' print str(DiffOperator) + ' test' print '==============' for i in xrange(result.nbComponents): - assert np.allclose(res_d[i][ind], ref_d[i][ind], rtol=errX) + print 'err = O(h**order) =', err[i % dom.dimension] + print np.max(np.abs(res_d[i][ind] - ref_d[i][ind])) + assert np.allclose(res_d[i][ind], ref_d[i][ind], + rtol=err[i % dom.dimension]) + #Op.finalize() def velocity_f(res, x, y, z, t): @@ -83,11 +95,97 @@ def grad_velo(res, x, y, z, t): return res -def test_Curl(): +def test_CurlFD(): + from parmepy.methods_keys import SpaceDiscretisation + from parmepy.methods import FD_C_4 from parmepy.operator.differential import Curl - callOp(Curl, vorticity_f) + method = {SpaceDiscretisation: FD_C_4} + callOp(Curl, vorticity_f, method=method, resolution=resol) + + +def test_CurlFD2(): + from parmepy.methods_keys import SpaceDiscretisation + from parmepy.methods import FD_C_2 + from parmepy.operator.differential import Curl + method = {SpaceDiscretisation: FD_C_2} + callOp(Curl, vorticity_f, method=method, order=2, resolution=resol) + + +def test_CurlFFT(): + from parmepy.methods_keys import SpaceDiscretisation + from parmepy.operator.differential import Curl + method = {SpaceDiscretisation: 'fftw'} + callOp(Curl, vorticity_f, method=method, order=6, resolution=resol) def test_Grad(): + from parmepy.methods_keys import SpaceDiscretisation + from parmepy.methods import FD_C_2 + method = {SpaceDiscretisation: FD_C_2} + from parmepy.operator.differential import Grad + callOp(Grad, grad_velo, op_dim=9, method=method, order=2, resolution=resol) + + +def test_Grad2(): + from parmepy.operator.differential import Grad + from parmepy.methods_keys import SpaceDiscretisation + from parmepy.methods import FD_C_4 + method = {SpaceDiscretisation: FD_C_4} + callOp(Grad, grad_velo, op_dim=9, method=method, resolution=resol) + + +def test_CurlFD_2(): + from parmepy.methods_keys import SpaceDiscretisation + from parmepy.methods import FD_C_4 + from parmepy.operator.differential import Curl + method = {SpaceDiscretisation: FD_C_4} + callOp(Curl, vorticity_f, method=method, dom=box2, resolution=resol2) + + +def test_CurlFD2_2(): + from parmepy.methods_keys import SpaceDiscretisation + from parmepy.methods import FD_C_2 + from parmepy.operator.differential import Curl + method = {SpaceDiscretisation: FD_C_2} + callOp(Curl, vorticity_f, method=method, order=2, dom=box2, + resolution=resol2) + + +def test_CurlFFT_2(): + from parmepy.methods_keys import SpaceDiscretisation + from parmepy.operator.differential import Curl + method = {SpaceDiscretisation: 'fftw'} + callOp(Curl, vorticity_f, method=method, order=6, dom=box2, + resolution=resol2) + + +def test_Grad_2(): + from parmepy.methods_keys import SpaceDiscretisation + from parmepy.methods import FD_C_2 + method = {SpaceDiscretisation: FD_C_2} + from parmepy.operator.differential import Grad + callOp(Grad, grad_velo, op_dim=9, method=method, order=2, dom=box2, + resolution=resol2) + + +def test_Grad2_2(): from parmepy.operator.differential import Grad - callOp(Grad, grad_velo, op_dim=9) + from parmepy.methods_keys import SpaceDiscretisation + from parmepy.methods import FD_C_4 + method = {SpaceDiscretisation: FD_C_4} + callOp(Grad, grad_velo, op_dim=9, method=method, dom=box2, + resolution=resol2) + + +# This may be useful to run mpi tests +if __name__ == "__main__": + test_CurlFD() + test_CurlFD2() + test_CurlFFT() + test_Grad() + test_Grad2() + test_CurlFD_2() + test_CurlFD2_2() + test_CurlFFT_2() + test_Grad_2() + test_Grad2_2() diff --git a/HySoP/hysop/operator/tests/test_particle_advection.py b/HySoP/hysop/operator/tests/test_particle_advection.py index 3a846186ecaf10626adfa03ef7a76b295b53500d..34773a79b5fa34fc8afba6849c1fb24fbe1981bf 100644 --- a/HySoP/hysop/operator/tests/test_particle_advection.py +++ b/HySoP/hysop/operator/tests/test_particle_advection.py @@ -64,7 +64,6 @@ def setup_list_3D(): def assertion(scal, advec): advec.discretize() advec.setUp() - scal_d = scal.discreteFields.values()[0] scal_d.data[0][...] = np.asarray(np.random.random(scal_d.data[0].shape), dtype=PARMES_REAL, order=ORDER) @@ -78,6 +77,7 @@ def assertion_vector2D(scal, advec): advec.discretize() advec.setUp() + print "aaaa ", advec.velocity.norm() scal_d = scal.discreteFields.values()[0] scal_d.data[0][...] = np.asarray(np.random.random(scal_d.data[0].shape), dtype=PARMES_REAL, order=ORDER) @@ -87,6 +87,8 @@ def assertion_vector2D(scal, advec): scal2_init = npw.copy(scal_d.data[1]) advec.apply(Simulation(tinit=0., tend=0.01, nbIter=1)) + print np.max(np.abs((scal1_init - scal_d.data[0]))) + print np.max(np.abs((scal2_init - scal_d.data[1]))) return np.allclose(scal1_init, scal_d.data[0]) and \ np.allclose(scal2_init, scal_d.data[1]) @@ -135,6 +137,7 @@ def assertion_list(scal, advec): def test_nullVelocity_2D(): """ """ + print "START ..." scal, velo = setup_2D() advec = Advection(velo, scal, @@ -207,3 +210,7 @@ def test_nullVelocity_vector_3D(): if __name__ == '__main__': test_nullVelocity_2D() + #test_nullVelocity_vector_2D() + test_nullVelocity_list_2D() + test_nullVelocity_3D() + test_nullVelocity_vector_3D() diff --git a/HySoP/hysop/problem/simulation.py b/HySoP/hysop/problem/simulation.py index 20247b3a14f9e20b4e0db8782c7965bae1de6950..309e8f175862be78a722d72b6523fcb16ef2abdd 100644 --- a/HySoP/hysop/problem/simulation.py +++ b/HySoP/hysop/problem/simulation.py @@ -46,7 +46,6 @@ class Simulation(object): print 'Warning : both nbIter and timeStep are given.\ timeStep is ignored' self.timeStep = (self.end - self.start) / self.nbIter - print "{0:24.16f}".format(self.timeStep) elif timeStep is not None: ## Simulation time step self.timeStep = timeStep