File: genshtmat.f90

package info (click to toggle)
elkcode 5.4.24-2
  • links: PTS, VCS
  • area: main
  • in suites: buster
  • size: 12,840 kB
  • sloc: f90: 48,415; fortran: 22,457; perl: 965; makefile: 384; sh: 369; python: 105; ansic: 67
file content (119 lines) | stat: -rw-r--r-- 3,626 bytes parent folder | download | duplicates (2)
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119

! Copyright (C) 2002-2005 J. K. Dewhurst, S. Sharma and C. Ambrosch-Draxl.
! This file is distributed under the terms of the GNU General Public License.
! See the file COPYING for license details.

!BOP
! !ROUTINE: genshtmat
! !INTERFACE:
subroutine genshtmat
! !USES:
use modmain
! !DESCRIPTION:
!   Generates the forward and backward spherical harmonic transformation (SHT)
!   matrices using the spherical covering set produced by the routine
!   {\tt sphcover}. These matrices are used to transform a function between its
!   $(l,m)$-expansion coefficients and its values at the $(\theta,\phi)$ points
!   on the sphere.
!
! !REVISION HISTORY:
!   Created April 2003 (JKD)
!EOP
!BOC
implicit none
! local variables
integer itp,info
real(8) v1(3),v2(3),r
! automatic arrays
integer ipiv(lmmaxo)
real(8) tp(2,lmmaxo),rlm(lmmaxo),work(lmmaxo)
complex(8) ylm(lmmaxo),zwork(lmmaxo)
!--------------------------------!
!     SHT matrices for lmaxo     !
!--------------------------------!
! allocate real SHT matrices
if (allocated(rbshto)) deallocate(rbshto)
allocate(rbshto(lmmaxo,lmmaxo))
if (allocated(rfshto)) deallocate(rfshto)
allocate(rfshto(lmmaxo,lmmaxo))
! allocate complex SHT matrices
if (allocated(zbshto)) deallocate(zbshto)
allocate(zbshto(lmmaxo,lmmaxo))
if (allocated(zfshto)) deallocate(zfshto)
allocate(zfshto(lmmaxo,lmmaxo))
! generate spherical covering set
call sphcover(lmmaxo,tp)
! rotate the spherical covering set if required
if (trotsht) then
  do itp=1,lmmaxo
    call sctovec(tp(:,itp),v1)
    call r3mv(rotsht,v1,v2)
    call sphcrd(v2,r,tp(:,itp))
  end do
end if
! generate real and complex spherical harmonics and set the backward SHT arrays
do itp=1,lmmaxo
  call genrlm(lmaxo,tp(:,itp),rlm)
  rbshto(itp,1:lmmaxo)=rlm(1:lmmaxo)
  call genylm(lmaxo,tp(:,itp),ylm)
  zbshto(itp,1:lmmaxo)=ylm(1:lmmaxo)
end do
! find the forward SHT arrays
! real
rfshto(:,:)=rbshto(:,:)
call dgetrf(lmmaxo,lmmaxo,rfshto,lmmaxo,ipiv,info)
if (info.ne.0) goto 10
call dgetri(lmmaxo,rfshto,lmmaxo,ipiv,work,lmmaxo,info)
if (info.ne.0) goto 10
! complex
zfshto(:,:)=zbshto(:,:)
call zgetrf(lmmaxo,lmmaxo,zfshto,lmmaxo,ipiv,info)
if (info.ne.0) goto 10
call zgetri(lmmaxo,zfshto,lmmaxo,ipiv,zwork,lmmaxo,info)
if (info.ne.0) goto 10
!--------------------------------!
!     SHT matrices for lmaxi     !
!--------------------------------!
! allocate real SHT matrices
if (allocated(rbshti)) deallocate(rbshti)
allocate(rbshti(lmmaxi,lmmaxi))
if (allocated(rfshti)) deallocate(rfshti)
allocate(rfshti(lmmaxi,lmmaxi))
! allocate complex SHT matrices
if (allocated(zbshti)) deallocate(zbshti)
allocate(zbshti(lmmaxi,lmmaxi))
if (allocated(zfshti)) deallocate(zfshti)
allocate(zfshti(lmmaxi,lmmaxi))
! generate spherical covering set for lmaxi
call sphcover(lmmaxi,tp)
! generate real and complex spherical harmonics and set the backward SHT arrays
do itp=1,lmmaxi
  call genrlm(lmaxi,tp(:,itp),rlm)
  rbshti(itp,1:lmmaxi)=rlm(1:lmmaxi)
  call genylm(lmaxi,tp(:,itp),ylm)
  zbshti(itp,1:lmmaxi)=ylm(1:lmmaxi)
end do
! find the forward SHT arrays
! real
rfshti(:,:)=rbshti(:,:)
call dgetrf(lmmaxi,lmmaxi,rfshti,lmmaxi,ipiv,info)
if (info.ne.0) goto 10
call dgetri(lmmaxi,rfshti,lmmaxi,ipiv,work,lmmaxi,info)
if (info.ne.0) goto 10
! complex
zfshti(:,:)=zbshti(:,:)
call zgetrf(lmmaxi,lmmaxi,zfshti,lmmaxi,ipiv,info)
if (info.ne.0) goto 10
call zgetri(lmmaxi,zfshti,lmmaxi,ipiv,zwork,lmmaxi,info)
if (info.ne.0) goto 10
return
10 continue
write(*,*)
write(*,'("Error(genshtmat): unable to find inverse spherical harmonic &
 &transform")')
write(*,'(" => improper spherical covering")')
write(*,*)
stop
end subroutine
!EOC