File: AngularF.c

package info (click to toggle)
openmx 3.7.6-1
  • links: PTS, VCS
  • area: main
  • in suites: jessie, jessie-kfreebsd, stretch
  • size: 325,856 kB
  • ctags: 3,575
  • sloc: ansic: 152,655; f90: 2,080; python: 876; makefile: 675; sh: 25; perl: 18
file content (213 lines) | stat: -rw-r--r-- 6,177 bytes parent folder | download | duplicates (4)
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
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
/**********************************************************************
  AngularF.c:

     AngularF.c is a subroutine to calculate the spherical solid
     harmonic function specified by "(l,m)" at (Q,P).

  Log of AngularF.c:

     22/Nov/2001  Released by T.Ozaki

***********************************************************************/

#include <stdio.h>
#include <math.h>
#include "openmx_common.h"

double AngularF(int l, int m, double Q, double P,
                int Use_switch,
                double siQ, double coQ,
                double siP, double coP)
{
  double Y;

  if (Use_switch==0){
    siQ = sin(Q); 
    coQ = cos(Q);
    siP = sin(P);
    coP = cos(P);
  }

  switch(l){

    /************
          s
    ************/

    case 0:

      /****************************************************
                     s,  Y = 0.50/sqrt(PI)
      ****************************************************/

      Y = 0.50/sqrt(PI);

    break;

    /************
          p
    ************/

    case 1: 
      switch(m){
        case 0:

          /****************************************************
                     x, 0.50*sqrt(3.0/PI)*sin(Q)*cos(P)
          ****************************************************/

          Y = 0.50*sqrt(3.0/PI)*siQ*coP;

        break;
        case 1:

          /****************************************************
                     y, 0.50*sqrt(3.0/PI)*sin(Q)*sin(P)
          ****************************************************/

          Y = 0.50*sqrt(3.0/PI)*siQ*siP;

        break;
        case 2:

          /****************************************************
                     z, 0.50*sqrt(3.0/PI)*cos(Q)
          ****************************************************/

          Y = 0.50*sqrt(3.0/PI)*coQ;
        break;
      }
    break;

    /************
          d
    ************/

    case 2: 
      switch(m){
        case 0:

          /****************************************************
             3zz-rr, sqrt(5.0/16.0/PI)*(3.0*cos(Q)*cos(Q)-1.0)
          ****************************************************/

          Y = 0.94617469575756*coQ*coQ - 0.31539156525252;
        break;
        case 1:

          /*******************************************************
                                   xx-yy,
                            2.0*sqrt(15.0/32.0/PI)*
                       sin(Q)*sin(Q)*cos(2.0*P)/sqrt(2.0)
          *******************************************************/

          Y = 0.54627421529604*siQ*siQ*(1.0 - 2.0*siP*siP);
        break;
        case 2:

          /*******************************************************
                                    xy,
                         2.0*sqrt(15.0/32.0/PI)*
                     sin(Q)*sin(Q)*sin(2.0*P)/sqrt(2.0)
          *******************************************************/

          Y = 1.09254843059208*siQ*siQ*siP*coP;
        break;
        case 3:

          /*******************************************************
                                    xz,
                           2.0*sqrt(15.0/8.0/PI)*
                        sin(Q)*cos(Q)*cos(P)/sqrt(2.0)
          *******************************************************/

          Y = 1.09254843059208*siQ*coQ*coP;
        break;
        case 4:

          /*******************************************************
                                    yz,
                           2.0*sqrt(15.0/8.0/PI)*
                       sin(Q)*cos(Q)*sin(P)/sqrt(2.0)
          *******************************************************/

          Y = 1.09254843059208*siQ*coQ*siP;
        break;
      }
    break;

    /************
          f
    ************/

    case 3:
      switch(m){
        case 0:

          /*******************************************************
                                 5z^2-3r^2,
                   sqrt(7)/4/sqrt(PI)*(5*cos(Q)^3-3*cos(Q))
          *******************************************************/

          Y = 0.373176332590116*(5.0*coQ*coQ*coQ - 3.0*coQ);
        break;
        case 1:

          /*******************************************************
                                 5xz^2-xr^2, 
               sqrt(42)/8/sqrt(PI)*cos(P)sin(Q)*(5*cos(Q)^2-1)
          *******************************************************/

          Y = 0.457045799464466*coP*siQ*(5.0*coQ*coQ - 1.0);
        break;
        case 2:

          /*******************************************************
                                 5yz^2-yr^2,
               sqrt(42)/8/sqrt(PI)*sin(P)sin(Q)*(5*cos(Q)^2-1)
          *******************************************************/

          Y = 0.457045799464466*siP*siQ*(5.0*coQ*coQ - 1.0);
        break;
        case 3:

          /*******************************************************
                                zx^2-zy^2,
                           sqrt(105)/4/sqrt(PI)*
                    sin(Q)^2cos(Q)*(cos(P)^2-sin(P)^2)
          *******************************************************/

          Y = 1.44530572132028*siQ*siQ*coQ*(coP*coP-siP*siP);
        break;
        case 4:

          /*******************************************************
                                    xyz,
              sqrt(105)/2/sqrt(PI)*sin(Q)^2cos(Q)*sin(P)*cos(P)
          *******************************************************/

          Y = 2.89061144264055*siQ*siQ*coQ*siP*coP;
        break;
        case 5:

          /*******************************************************
                                 x^3-3*xy^2,
              sqrt(70)/8/sqrt(PI)*sin(Q)^3*(4*cos(P)^3-3*cos(P))
          *******************************************************/

          Y = 0.590043589926644*siQ*siQ*siQ*(4.0*coP*coP*coP - 3.0*coP);
        break;
        case 6:

          /*******************************************************
                                  3yx^2-y^3,
             sqrt(70)/8/sqrt(PI)*sin(Q)^3*(3*sin(P)-4*sin(P)^3)
          *******************************************************/

          Y = 0.590043589926644*siQ*siQ*siQ*(3.0*siP - 4.0*siP*siP*siP);
        break;
      }
    break;
  }
  return Y;
}