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-2022 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 .

SOURCE

17 #if defined HAVE_CONFIG_H
18 #include "config.h"
19 #endif
20 
21 #include "abi_common.h"
22 
23 MODULE m_lebedev
24 
25  use defs_basis
26  use m_abicore
27  use m_errors
28 
29  use m_fstrings,     only : sjoin, itoa
30 
31  implicit none
32 
33  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.

SOURCE

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

m_lebedev/gen_oh [ Functions ]

[ Top ] [ m_lebedev ] [ Functions ]

NAME

  gen_oh

FUNCTION

INPUTS

OUTPUT

SOURCE

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

m_lebedev/lebedev_free [ Functions ]

[ Top ] [ m_lebedev ] [ Functions ]

NAME

  lebedev_free

FUNCTION

  Free an instance of lebedev_t

SOURCE

215 subroutine lebedev_free(lgrid)
216 
217 !Arguments ------------------------------------
218  class(lebedev_t),intent(inout) :: lgrid
219 
220 ! *********************************************************************
221 
222  lgrid%npts=0
223  ABI_SFREE(lgrid%versors)
224  ABI_SFREE(lgrid%weights)
225 
226 end subroutine lebedev_free

m_lebedev/lebedev_from_index [ Functions ]

[ Top ] [ m_lebedev ] [ Functions ]

NAME

  lebedev_from_index

FUNCTION

  Create a lebedev grid from the sequential index.
  Useful if one has to loop over all the grids.

INPUTS

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

OUTPUT

  Lgrid<lebedev_t>=The grid fully initialized.

SOURCE

178 subroutine lebedev_from_index(new, seq_idx)
179 
180 !Arguments ------------------------------------
181 !scalars
182  class(lebedev_t),intent(out) :: new
183  integer,intent(in) :: seq_idx
184 
185 !Local variables-------------------------------
186  integer :: ii
187  real(dp),allocatable :: xx(:),yy(:),zz(:)
188 ! *********************************************************************
189 
190  call build_lebedev_grid(seq_idx, new%npts, xx, yy, zz, new%weights)
191 
192  ABI_MALLOC(new%versors, (3,new%npts))
193  do ii=1,new%npts
194    new%versors(:,ii) = [xx(ii), yy(ii), zz(ii)]
195  end do
196 
197  ABI_FREE(xx)
198  ABI_FREE(yy)
199  ABI_FREE(zz)
200 
201 end subroutine lebedev_from_index

m_lebedev/lebedev_from_npts [ Functions ]

[ Top ] [ m_lebedev ] [ Functions ]

NAME

  lebedev_from_npts

FUNCTION

  Create a lebedev grid given the number of points. npts
  Exit status:
    0 on success
   -1 if npts is greater than the max number of points supported.
   N > 0 is the first grid with size greater that npts.

INPUTS

SOURCE

134 subroutine lebedev_from_npts(new, npts, ierr)
135 
136 !Arguments ------------------------------------
137 !scalars
138  class(lebedev_t),intent(out) :: new
139  integer,intent(in) :: npts
140  integer,intent(out) :: ierr
141 
142 !Local variables-------------------------------
143  integer :: seq_idx
144 ! *********************************************************************
145 
146  ierr = -1
147  do seq_idx=1, lebedev_ngrids
148    if (lebedev_npts(seq_idx) == npts) then
149      call new%from_index(seq_idx)
150      ierr = 0
151      return
152    else if (lebedev_npts(seq_idx) > npts .and. ierr == -1) then
153      ierr = lebedev_npts(seq_idx)
154    end if
155  end do
156 
157 end subroutine lebedev_from_npts

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

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