TABLE OF CONTENTS


ABINIT/m_lebedev [ Modules ]

[ Top ] [ Modules ]

NAME

  m_lebedev

FUNCTION

  This module contains routines for performing integrations
  on the sphere using lebedev-laikov angular grids.

COPYRIGHT

 Copyright (C) 2008-2018 ABINIT group (MG)
 This file is distributed under the terms of the
 GNU General Public License, see ~abinit/COPYING
 or http://www.gnu.org/copyleft/gpl.txt .

PARENTS

CHILDREN

SOURCE

22 #if defined HAVE_CONFIG_H
23 #include "config.h"
24 #endif
25 
26 #include "abi_common.h"
27 
28 MODULE m_lebedev
29 
30  use defs_basis
31  use m_abicore
32  use m_errors
33 
34  use m_fstrings,     only : sjoin, itoa
35 
36  implicit none
37 
38  private

m_lebedev/build_lebedev_grid [ Functions ]

[ Top ] [ m_lebedev ] [ Functions ]

NAME

  build_lebedev_grid

FUNCTION

  Helper subroutine returning the knots and weights of the angular grid
  from its sequential index

INPUTS

  seq_idx=sequential index (must be in [1,32]

OUTPUT

  npts=Number of points
  xx(npts),yy(npts),zz(npts)=The Cartesian coordinates of the knots.
  ww(npts)=The weights.

PARENTS

      m_lebedev

CHILDREN

      gen_oh

SOURCE

251 subroutine build_lebedev_grid(seq_idx, npts, xx, yy, zz, ww)
252 
253 
254 !This section has been created automatically by the script Abilint (TD).
255 !Do not modify the following lines by hand.
256 #undef ABI_FUNC
257 #define ABI_FUNC 'build_lebedev_grid'
258 !End of the abilint section
259 
260  implicit none
261 
262 !This section has been created automatically by the script Abilint (TD).
263 !Do not modify the following lines by hand.
264 #undef ABI_FUNC
265 #define ABI_FUNC 'build_lebedev_grid'
266 !End of the abilint section
267 
268 !Arguments ------------------------------------
269 !scalars
270  integer,intent(in) :: seq_idx
271  integer,intent(out) :: npts
272  real(dp),allocatable,intent(out) :: xx(:),yy(:),zz(:),ww(:)
273 
274 ! *********************************************************************
275 
276  if (seq_idx < 1 .or. seq_idx > lebedev_ngrids) then
277    MSG_ERROR(sjoin("seq_idx ",itoa(seq_idx), "out of range"))
278  end if
279 
280  npts = lebedev_npts(seq_idx)
281  ABI_MALLOC(xx, (npts))
282  ABI_MALLOC(yy, (npts))
283  ABI_MALLOC(zz, (npts))
284  ABI_MALLOC(ww, (npts))
285 
286  SELECT CASE (seq_idx)
287  CASE (1)
288    call LD0006(xx,yy,zz,ww,npts)
289  CASE (2)
290    call LD0014(xx,yy,zz,ww,npts)
291  CASE (3)
292    call LD0026(xx,yy,zz,ww,npts)
293  CASE (4)
294    call LD0038(xx,yy,zz,ww,npts)
295  CASE (5)
296    call LD0050(xx,yy,zz,ww,npts)
297  CASE (6)
298    call LD0074(xx,yy,zz,ww,npts)
299  CASE (7)
300    call LD0086(xx,yy,zz,ww,npts)
301  CASE (8)
302    call LD0110(xx,yy,zz,ww,npts)
303  CASE (9)
304    call LD0146(xx,yy,zz,ww,npts)
305  CASE (10)
306    call LD0170(xx,yy,zz,ww,npts)
307  CASE (11)
308    call LD0194(xx,yy,zz,ww,npts)
309  CASE (12)
310    call LD0230(xx,yy,zz,ww,npts)
311  CASE (13)
312    call LD0266(xx,yy,zz,ww,npts)
313  CASE (14)
314    call LD0302(xx,yy,zz,ww,npts)
315  CASE (15)
316    call LD0350(xx,yy,zz,ww,npts)
317  CASE (16)
318    call LD0434(xx,yy,zz,ww,npts)
319  CASE (17)
320    call LD0590(xx,yy,zz,ww,npts)
321  CASE (18)
322    call LD0770(xx,yy,zz,ww,npts)
323  CASE (19)
324    call LD0974(xx,yy,zz,ww,npts)
325  CASE (20)
326    call LD1202(xx,yy,zz,ww,npts)
327  CASE (21)
328    call LD1454(xx,yy,zz,ww,npts)
329  CASE (22)
330    call LD1730(xx,yy,zz,ww,npts)
331  CASE (23)
332    call LD2030(xx,yy,zz,ww,npts)
333  CASE (24)
334    call LD2354(xx,yy,zz,ww,npts)
335  CASE (25)
336    call LD2702(xx,yy,zz,ww,npts)
337  CASE (26)
338    call LD3074(xx,yy,zz,ww,npts)
339  CASE (27)
340    call LD3470(xx,yy,zz,ww,npts)
341  CASE (28)
342    call LD3890(xx,yy,zz,ww,npts)
343  CASE (29)
344    call LD4334(xx,yy,zz,ww,npts)
345  CASE (30)
346    call LD4802(xx,yy,zz,ww,npts)
347  CASE (31)
348    call LD5294(xx,yy,zz,ww,npts)
349  CASE (32)
350    call LD5810(xx,yy,zz,ww,npts)
351  CASE DEFAULT
352    MSG_ERROR(sjoin("seq_idx ", itoa(seq_idx), "out of range"))
353  END SELECT
354 
355 end subroutine build_lebedev_grid

m_lebedev/gen_oh [ Functions ]

[ Top ] [ m_lebedev ] [ Functions ]

NAME

  gen_oh

FUNCTION

INPUTS

OUTPUT

PARENTS

      m_lebedev

CHILDREN

      gen_oh

SOURCE

378 subroutine gen_oh(code, num, x, y, z, w, a, b, v)
379 
380 
381 !This section has been created automatically by the script Abilint (TD).
382 !Do not modify the following lines by hand.
383 #undef ABI_FUNC
384 #define ABI_FUNC 'gen_oh'
385 !End of the abilint section
386 
387        implicit logical(a-z)
388        double precision x(*),y(*),z(*),w(*)
389        double precision a,b,v
390        integer code
391        integer num
392        double precision c
393 !chvd
394 !chvd   This subroutine is part of a set of subroutines that generate
395 !chvd   Lebedev grids [1-6] for integration on a sphere. The original
396 !chvd   C-code [1] was kindly provided by Dr. Dmitri N. Laikov and
397 !chvd   translated into fortran by Dr. Christoph van Wuellen.
398 !chvd   This subroutine was translated from C to fortran77 by hand.
399 !chvd
400 !chvd   Users of this code are asked to include reference [1] in their
401 !chvd   publications, and in the user- and programmers-manuals
402 !chvd   describing their codes.
403 !chvd
404 !chvd   This code was distributed through CCL (http://www.ccl.net/).
405 !chvd
406 !chvd   [1] V.I. Lebedev, and D.N. Laikov
407 !chvd       "A quadrature formula for the sphere of the 131st
408 !chvd        algebraic order of accuracy"
409 !chvd       Doklady Mathematics, Vol. 59, No. 3, 1999, pp. 477-481. [[cite:Lebedev1999]]
410 !chvd
411 !chvd   [2] V.I. Lebedev
412 !chvd       "A quadrature formula for the sphere of 59th algebraic
413 !chvd        order of accuracy"
414 !chvd       Russian Acad. Sci. Dokl. Math., Vol. 50, 1995, pp. 283-286. [[cite:Lebedev1995]]
415 !chvd
416 !chvd   [3] V.I. Lebedev, and A.L. Skorokhodov
417 !chvd       "Quadrature formulas of orders 41, 47, and 53 for the sphere"
418 !chvd       Russian Acad. Sci. Dokl. Math., Vol. 45, 1992, pp. 587-592. [[cite:Lebedev1992]]
419 !chvd
420 !chvd   [4] V.I. Lebedev
421 !chvd       "Spherical quadrature formulas exact to orders 25-29"
422 !chvd       Siberian Mathematical Journal, Vol. 18, 1977, pp. 99-107. [[cite:Lebedev1977]]
423 !chvd
424 !chvd   [5] V.I. Lebedev
425 !chvd       "Quadratures on a sphere"
426 !chvd       Computational Mathematics and Mathematical Physics, Vol. 16, 1976, pp. 10-24. [[cite:Lebedev1976]]
427 !chvd
428 !chvd   [6] V.I. Lebedev
429 !chvd       "Values of the nodes and weights of ninth to seventeenth
430 !chvd        order Gauss-Markov quadrature formulae invariant under the
431 !chvd        octahedron group with inversion"
432 !chvd       Computational Mathematics and Mathematical Physics, Vol. 15, 1975, pp. 44-51. [[cite:Lebedev1975]]
433 !chvd
434 !cvw
435 !cvw    Given a point on a sphere (specified by a and b), generate all
436 !cvw    the equivalent points under Oh symmetry, making grid points with
437 !cvw    weight v.
438 !cvw    The variable num is increased by the number of different points
439 !cvw    generated.
440 !cvw
441 !cvw    Depending on code, there are 6...48 different but equivalent
442 !cvw    points.
443 !cvw
444 !cvw    code=1:   (0,0,1) etc                                (  6 points)
445 !cvw    code=2:   (0,a,a) etc, a=1/sqrt(2)                   ( 12 points)
446 !cvw    code=3:   (a,a,a) etc, a=1/sqrt(3)                   (  8 points)
447 !cvw    code=4:   (a,a,b) etc, b=sqrt(1-2 a^2)               ( 24 points)
448 !cvw    code=5:   (a,b,0) etc, b=sqrt(1-a^2), a input        ( 24 points)
449 !cvw    code=6:   (a,b,c) etc, c=sqrt(1-a^2-b^2), a/b input  ( 48 points)
450 !cvw
451        select case(code)
452        case(1)
453        a=1.0d0
454        x(1) =  a
455        y(1) =  0.0d0
456        z(1) =  0.0d0
457        w(1) =  v
458        x(2) = -a
459        y(2) =  0.0d0
460        z(2) =  0.0d0
461        w(2) =  v
462        x(3) =  0.0d0
463        y(3) =  a
464        z(3) =  0.0d0
465        w(3) =  v
466        x(4) =  0.0d0
467        y(4) = -a
468        z(4) =  0.0d0
469        w(4) =  v
470        x(5) =  0.0d0
471        y(5) =  0.0d0
472        z(5) =  a
473        w(5) =  v
474        x(6) =  0.0d0
475        y(6) =  0.0d0
476        z(6) = -a
477        w(6) =  v
478        num=num+6
479 !cvw
480        case(2)
481        a=sqrt(0.5d0)
482        x( 1) =  0d0
483        y( 1) =  a
484        z( 1) =  a
485        w( 1) =  v
486        x( 2) =  0d0
487        y( 2) = -a
488        z( 2) =  a
489        w( 2) =  v
490        x( 3) =  0d0
491        y( 3) =  a
492        z( 3) = -a
493        w( 3) =  v
494        x( 4) =  0d0
495        y( 4) = -a
496        z( 4) = -a
497        w( 4) =  v
498        x( 5) =  a
499        y( 5) =  0d0
500        z( 5) =  a
501        w( 5) =  v
502        x( 6) = -a
503        y( 6) =  0d0
504        z( 6) =  a
505        w( 6) =  v
506        x( 7) =  a
507        y( 7) =  0d0
508        z( 7) = -a
509        w( 7) =  v
510        x( 8) = -a
511        y( 8) =  0d0
512        z( 8) = -a
513        w( 8) =  v
514        x( 9) =  a
515        y( 9) =  a
516        z( 9) =  0d0
517        w( 9) =  v
518        x(10) = -a
519        y(10) =  a
520        z(10) =  0d0
521        w(10) =  v
522        x(11) =  a
523        y(11) = -a
524        z(11) =  0d0
525        w(11) =  v
526        x(12) = -a
527        y(12) = -a
528        z(12) =  0d0
529        w(12) =  v
530        num=num+12
531 !cvw
532        case(3)
533        a = sqrt(1d0/3d0)
534        x(1) =  a
535        y(1) =  a
536        z(1) =  a
537        w(1) =  v
538        x(2) = -a
539        y(2) =  a
540        z(2) =  a
541        w(2) =  v
542        x(3) =  a
543        y(3) = -a
544        z(3) =  a
545        w(3) =  v
546        x(4) = -a
547        y(4) = -a
548        z(4) =  a
549        w(4) =  v
550        x(5) =  a
551        y(5) =  a
552        z(5) = -a
553        w(5) =  v
554        x(6) = -a
555        y(6) =  a
556        z(6) = -a
557        w(6) =  v
558        x(7) =  a
559        y(7) = -a
560        z(7) = -a
561        w(7) =  v
562        x(8) = -a
563        y(8) = -a
564        z(8) = -a
565        w(8) =  v
566        num=num+8
567 !cvw
568        case(4)
569        b = sqrt(1d0 - 2d0*a*a)
570        x( 1) =  a
571        y( 1) =  a
572        z( 1) =  b
573        w( 1) =  v
574        x( 2) = -a
575        y( 2) =  a
576        z( 2) =  b
577        w( 2) =  v
578        x( 3) =  a
579        y( 3) = -a
580        z( 3) =  b
581        w( 3) =  v
582        x( 4) = -a
583        y( 4) = -a
584        z( 4) =  b
585        w( 4) =  v
586        x( 5) =  a
587        y( 5) =  a
588        z( 5) = -b
589        w( 5) =  v
590        x( 6) = -a
591        y( 6) =  a
592        z( 6) = -b
593        w( 6) =  v
594        x( 7) =  a
595        y( 7) = -a
596        z( 7) = -b
597        w( 7) =  v
598        x( 8) = -a
599        y( 8) = -a
600        z( 8) = -b
601        w( 8) =  v
602        x( 9) =  a
603        y( 9) =  b
604        z( 9) =  a
605        w( 9) =  v
606        x(10) = -a
607        y(10) =  b
608        z(10) =  a
609        w(10) =  v
610        x(11) =  a
611        y(11) = -b
612        z(11) =  a
613        w(11) =  v
614        x(12) = -a
615        y(12) = -b
616        z(12) =  a
617        w(12) =  v
618        x(13) =  a
619        y(13) =  b
620        z(13) = -a
621        w(13) =  v
622        x(14) = -a
623        y(14) =  b
624        z(14) = -a
625        w(14) =  v
626        x(15) =  a
627        y(15) = -b
628        z(15) = -a
629        w(15) =  v
630        x(16) = -a
631        y(16) = -b
632        z(16) = -a
633        w(16) =  v
634        x(17) =  b
635        y(17) =  a
636        z(17) =  a
637        w(17) =  v
638        x(18) = -b
639        y(18) =  a
640        z(18) =  a
641        w(18) =  v
642        x(19) =  b
643        y(19) = -a
644        z(19) =  a
645        w(19) =  v
646        x(20) = -b
647        y(20) = -a
648        z(20) =  a
649        w(20) =  v
650        x(21) =  b
651        y(21) =  a
652        z(21) = -a
653        w(21) =  v
654        x(22) = -b
655        y(22) =  a
656        z(22) = -a
657        w(22) =  v
658        x(23) =  b
659        y(23) = -a
660        z(23) = -a
661        w(23) =  v
662        x(24) = -b
663        y(24) = -a
664        z(24) = -a
665        w(24) =  v
666        num=num+24
667 !cvw
668        case(5)
669        b=sqrt(1d0-a*a)
670        x( 1) =  a
671        y( 1) =  b
672        z( 1) =  0d0
673        w( 1) =  v
674        x( 2) = -a
675        y( 2) =  b
676        z( 2) =  0d0
677        w( 2) =  v
678        x( 3) =  a
679        y( 3) = -b
680        z( 3) =  0d0
681        w( 3) =  v
682        x( 4) = -a
683        y( 4) = -b
684        z( 4) =  0d0
685        w( 4) =  v
686        x( 5) =  b
687        y( 5) =  a
688        z( 5) =  0d0
689        w( 5) =  v
690        x( 6) = -b
691        y( 6) =  a
692        z( 6) =  0d0
693        w( 6) =  v
694        x( 7) =  b
695        y( 7) = -a
696        z( 7) =  0d0
697        w( 7) =  v
698        x( 8) = -b
699        y( 8) = -a
700        z( 8) =  0d0
701        w( 8) =  v
702        x( 9) =  a
703        y( 9) =  0d0
704        z( 9) =  b
705        w( 9) =  v
706        x(10) = -a
707        y(10) =  0d0
708        z(10) =  b
709        w(10) =  v
710        x(11) =  a
711        y(11) =  0d0
712        z(11) = -b
713        w(11) =  v
714        x(12) = -a
715        y(12) =  0d0
716        z(12) = -b
717        w(12) =  v
718        x(13) =  b
719        y(13) =  0d0
720        z(13) =  a
721        w(13) =  v
722        x(14) = -b
723        y(14) =  0d0
724        z(14) =  a
725        w(14) =  v
726        x(15) =  b
727        y(15) =  0d0
728        z(15) = -a
729        w(15) =  v
730        x(16) = -b
731        y(16) =  0d0
732        z(16) = -a
733        w(16) =  v
734        x(17) =  0d0
735        y(17) =  a
736        z(17) =  b
737        w(17) =  v
738        x(18) =  0d0
739        y(18) = -a
740        z(18) =  b
741        w(18) =  v
742        x(19) =  0d0
743        y(19) =  a
744        z(19) = -b
745        w(19) =  v
746        x(20) =  0d0
747        y(20) = -a
748        z(20) = -b
749        w(20) =  v
750        x(21) =  0d0
751        y(21) =  b
752        z(21) =  a
753        w(21) =  v
754        x(22) =  0d0
755        y(22) = -b
756        z(22) =  a
757        w(22) =  v
758        x(23) =  0d0
759        y(23) =  b
760        z(23) = -a
761        w(23) =  v
762        x(24) =  0d0
763        y(24) = -b
764        z(24) = -a
765        w(24) =  v
766        num=num+24
767 !cvw
768        case(6)
769        c=sqrt(1d0 - a*a - b*b)
770        x( 1) =  a
771        y( 1) =  b
772        z( 1) =  c
773        w( 1) =  v
774        x( 2) = -a
775        y( 2) =  b
776        z( 2) =  c
777        w( 2) =  v
778        x( 3) =  a
779        y( 3) = -b
780        z( 3) =  c
781        w( 3) =  v
782        x( 4) = -a
783        y( 4) = -b
784        z( 4) =  c
785        w( 4) =  v
786        x( 5) =  a
787        y( 5) =  b
788        z( 5) = -c
789        w( 5) =  v
790        x( 6) = -a
791        y( 6) =  b
792        z( 6) = -c
793        w( 6) =  v
794        x( 7) =  a
795        y( 7) = -b
796        z( 7) = -c
797        w( 7) =  v
798        x( 8) = -a
799        y( 8) = -b
800        z( 8) = -c
801        w( 8) =  v
802        x( 9) =  a
803        y( 9) =  c
804        z( 9) =  b
805        w( 9) =  v
806        x(10) = -a
807        y(10) =  c
808        z(10) =  b
809        w(10) =  v
810        x(11) =  a
811        y(11) = -c
812        z(11) =  b
813        w(11) =  v
814        x(12) = -a
815        y(12) = -c
816        z(12) =  b
817        w(12) =  v
818        x(13) =  a
819        y(13) =  c
820        z(13) = -b
821        w(13) =  v
822        x(14) = -a
823        y(14) =  c
824        z(14) = -b
825        w(14) =  v
826        x(15) =  a
827        y(15) = -c
828        z(15) = -b
829        w(15) =  v
830        x(16) = -a
831        y(16) = -c
832        z(16) = -b
833        w(16) =  v
834        x(17) =  b
835        y(17) =  a
836        z(17) =  c
837        w(17) =  v
838        x(18) = -b
839        y(18) =  a
840        z(18) =  c
841        w(18) =  v
842        x(19) =  b
843        y(19) = -a
844        z(19) =  c
845        w(19) =  v
846        x(20) = -b
847        y(20) = -a
848        z(20) =  c
849        w(20) =  v
850        x(21) =  b
851        y(21) =  a
852        z(21) = -c
853        w(21) =  v
854        x(22) = -b
855        y(22) =  a
856        z(22) = -c
857        w(22) =  v
858        x(23) =  b
859        y(23) = -a
860        z(23) = -c
861        w(23) =  v
862        x(24) = -b
863        y(24) = -a
864        z(24) = -c
865        w(24) =  v
866        x(25) =  b
867        y(25) =  c
868        z(25) =  a
869        w(25) =  v
870        x(26) = -b
871        y(26) =  c
872        z(26) =  a
873        w(26) =  v
874        x(27) =  b
875        y(27) = -c
876        z(27) =  a
877        w(27) =  v
878        x(28) = -b
879        y(28) = -c
880        z(28) =  a
881        w(28) =  v
882        x(29) =  b
883        y(29) =  c
884        z(29) = -a
885        w(29) =  v
886        x(30) = -b
887        y(30) =  c
888        z(30) = -a
889        w(30) =  v
890        x(31) =  b
891        y(31) = -c
892        z(31) = -a
893        w(31) =  v
894        x(32) = -b
895        y(32) = -c
896        z(32) = -a
897        w(32) =  v
898        x(33) =  c
899        y(33) =  a
900        z(33) =  b
901        w(33) =  v
902        x(34) = -c
903        y(34) =  a
904        z(34) =  b
905        w(34) =  v
906        x(35) =  c
907        y(35) = -a
908        z(35) =  b
909        w(35) =  v
910        x(36) = -c
911        y(36) = -a
912        z(36) =  b
913        w(36) =  v
914        x(37) =  c
915        y(37) =  a
916        z(37) = -b
917        w(37) =  v
918        x(38) = -c
919        y(38) =  a
920        z(38) = -b
921        w(38) =  v
922        x(39) =  c
923        y(39) = -a
924        z(39) = -b
925        w(39) =  v
926        x(40) = -c
927        y(40) = -a
928        z(40) = -b
929        w(40) =  v
930        x(41) =  c
931        y(41) =  b
932        z(41) =  a
933        w(41) =  v
934        x(42) = -c
935        y(42) =  b
936        z(42) =  a
937        w(42) =  v
938        x(43) =  c
939        y(43) = -b
940        z(43) =  a
941        w(43) =  v
942        x(44) = -c
943        y(44) = -b
944        z(44) =  a
945        w(44) =  v
946        x(45) =  c
947        y(45) =  b
948        z(45) = -a
949        w(45) =  v
950        x(46) = -c
951        y(46) =  b
952        z(46) = -a
953        w(46) =  v
954        x(47) =  c
955        y(47) = -b
956        z(47) = -a
957        w(47) =  v
958        x(48) = -c
959        y(48) = -b
960        z(48) = -a
961        w(48) =  v
962        num=num+48
963        case default
964        MSG_ERROR('Gen_Oh: Invalid Code')
965        stop
966        end select
967  end subroutine gen_oh

m_lebedev/lebedev_free [ Functions ]

[ Top ] [ m_lebedev ] [ Functions ]

NAME

  lebedev_free

FUNCTION

  Free an instance of lebedev_t

PARENTS

      m_ifc

CHILDREN

      gen_oh

SOURCE

191 subroutine lebedev_free(lgrid)
192 
193 
194 !This section has been created automatically by the script Abilint (TD).
195 !Do not modify the following lines by hand.
196 #undef ABI_FUNC
197 #define ABI_FUNC 'lebedev_free'
198 !End of the abilint section
199 
200  implicit none
201 
202 !This section has been created automatically by the script Abilint (TD).
203 !Do not modify the following lines by hand.
204 #undef ABI_FUNC
205 #define ABI_FUNC 'lebedev_free'
206 !End of the abilint section
207 
208 !Arguments ------------------------------------
209 !scalars
210  type(lebedev_t),intent(inout) :: lgrid
211 
212 ! *********************************************************************
213 
214  lgrid%npts=0
215  if (allocated(lgrid%versors)) then
216    ABI_FREE(lgrid%versors)
217  end if
218  if (allocated(lgrid%weights)) then
219    ABI_FREE(lgrid%weights)
220  end if
221 
222 end subroutine lebedev_free

m_lebedev/lebedev_new [ Functions ]

[ Top ] [ m_lebedev ] [ Functions ]

NAME

  lebedev_new

FUNCTION

  Create a lebedev grid.

INPUTS

  seq_idx=Sequential index comprised between 1 and 32 defining the order of the mesh.

OUTPUT

  Lgrid<lebedev_t>=The grid fully initialized.

PARENTS

      m_lebedev

CHILDREN

      gen_oh

SOURCE

134 type(lebedev_t) function lebedev_new(seq_idx) result(new)
135 
136 
137 !This section has been created automatically by the script Abilint (TD).
138 !Do not modify the following lines by hand.
139 #undef ABI_FUNC
140 #define ABI_FUNC 'lebedev_new'
141 !End of the abilint section
142 
143  implicit none
144 
145 !This section has been created automatically by the script Abilint (TD).
146 !Do not modify the following lines by hand.
147 #undef ABI_FUNC
148 #define ABI_FUNC 'lebedev_new'
149 !End of the abilint section
150 
151 !Arguments ------------------------------------
152 !scalars
153  integer,intent(in) :: seq_idx
154 
155 !Local variables-------------------------------
156  integer :: ii
157  real(dp),allocatable :: xx(:),yy(:),zz(:)
158 ! *********************************************************************
159 
160  call build_lebedev_grid(seq_idx,new%npts,xx,yy,zz,new%weights)
161 
162  ABI_MALLOC(new%versors,(3,new%npts))
163  do ii=1,new%npts
164    new%versors(:,ii) = [xx(ii), yy(ii), zz(ii)]
165  end do
166 
167  ABI_FREE(xx)
168  ABI_FREE(yy)
169  ABI_FREE(zz)
170 
171 end function lebedev_new

m_lebedev/lebedev_t [ Types ]

[ Top ] [ m_lebedev ] [ Types ]

NAME

 lebedev_t

FUNCTION

 Structure storing the knots and the weights of the lebedev-laikov grid.

SOURCE

 52  type,public :: lebedev_t
 53 
 54    integer :: npts
 55    ! Number of points in the grid
 56 
 57    real(dp),allocatable :: versors(:,:)
 58    ! versors(3, npts)
 59    ! Points on the sphere.
 60 
 61    real(dp),allocatable :: weights(:)
 62    ! weights(npts)
 63    ! weights for Spherical integration.
 64 
 65  end type lebedev_t
 66 
 67  public :: lebedev_new         ! Create a lebedev grid.
 68  public :: lebedev_free        ! Free memory
 69 
 70 !----------------------------------------------------------------------
 71 
 72  integer,public,parameter :: lebedev_ngrids=32
 73  ! The number of grids available.
 74 
 75  ! The number of points in each grid.
 76  integer,public,parameter :: lebedev_npts(lebedev_ngrids)=(/ &
 77 0006,&
 78 0014,&
 79 0026,&
 80 0038,&
 81 0050,&
 82 0074,&
 83 0086,&
 84 0110,&
 85 0146,&
 86 0170,&
 87 0194,&
 88 0230,&
 89 0266,&
 90 0302,&
 91 0350,&
 92 0434,&
 93 0590,&
 94 0770,&
 95 0974,&
 96 1202,&
 97 1454,&
 98 1730,&
 99 2030,&
100 2354,&
101 2702,&
102 3074,&
103 3470,&
104 3890,&
105 4334,&
106 4802,&
107 5294,&
108 5810/)
109 
110 contains  !===========================================================