Coverage for gpaw/test/poisson/test_screened_poisson.py: 100%

1 

2from gpaw.utilities.tools import coordinates 

3from gpaw.utilities.gauss import Gaussian 

4from gpaw.grid_descriptor import GridDescriptor 

5from gpaw.helmholtz import HelmholtzSolver, ScreenedPoissonGaussian 

6 

7# Initialize classes 

8 

9 

10def test_poisson_screened_poisson(): 

11 a = 20 # Size of cell 

12 inv_width = 31 # inverse width of the gaussian 

13 N = 48 # Number of grid points 

14 coupling = -0.4 # dampening 

15 Nc = (N, N, N) # Number of grid points along each axis 

16 gd = GridDescriptor(Nc, (a, a, a), 0) # Grid-descriptor object 

17 solver = HelmholtzSolver(k2=coupling, nn=3) # Numerical poisson solver 

18 # solver = PoissonSolver(nn=3) # Numerical poisson solver 

19 # solver = HelmholtzSolver(0.16) # Numerical poisson solver 

20 solver.set_grid_descriptor(gd) 

21 # Matrix with the square of the radial coordinate 

22 xyz, r2 = coordinates(gd) 

23 gauss = Gaussian(gd, a=inv_width) # An instance of Gaussian 

24 test_screened_poisson = ScreenedPoissonGaussian(gd, a=inv_width) 

25 

26 # Check if Gaussian potentials are made correctly 

27 

28 # Array for storing the potential 

29 pot = gd.zeros(dtype=float, global_array=False) 

30 solver.load_gauss() 

31 vg = test_screened_poisson.get_phi(-coupling) # esp. for dampening 

32 # Get analytic functions 

33 ng = gauss.get_gauss(0) 

34 # vg = solver.phi_gauss 

35 # Solve potential numerically 

36 solver.solve(pot, ng, charge=None, zero_initial_phi=True) 

37 # Determine residual 

38 # residual = norm(pot - vg) 

39 residual = gd.integrate((pot - vg)**2)**0.5 

40 

41 # print result 

42 print('residual %s' % ( 

43 residual)) 

44 assert residual < 0.003 

45 

46 # mpirun -np 2 python gauss_func.py --gpaw-parallel --gpaw-debug