@@ -71,160 +71,205 @@ def run(nlambda_list, ndigit_list):
71
71
ndigit_str = ' ' .join (map (str , ndigit_list ))
72
72
print (
73
73
f"""\
74
- ! Precomputed data of IR for a parameter matrix of
75
- ! nlambda = { nlambda_str } and ndigit = { ndigit_str } .
76
- ! This file was automatically generated:
77
- ! python3 mk_preset.py --nlambda { nlambda_str } --ndigit { ndigit_str } > sparse_ir_preset.f90
78
- ! Do not edit this file manually.
79
- !
80
- module sparse_ir_preset
81
- use sparse_ir
82
- implicit none
83
- """ )
74
+ !! Precomputed data of IR for a parameter matrix of
75
+ !! nlambda = { nlambda_str } and ndigit = { ndigit_str } .
76
+ !! This file was automatically generated:
77
+ !! python3 mk_preset.py --nlambda { nlambda_str } --ndigit { ndigit_str } > sparse_ir_preset.f90
78
+ !! Do not edit this file manually.
79
+ !!
80
+ !-----------------------------------------------------------------------
81
+ MODULE sparse_ir_preset
82
+ !-----------------------------------------------------------------------
83
+ USE sparse_ir
84
+ !
85
+ IMPLICIT NONE
86
+ !
87
+ PRIVATE
88
+ !
89
+ INTEGER, PARAMETER :: DP = KIND(0d0)
90
+ !
91
+ PUBLIC :: mk_ir_preset
92
+ !
93
+ """ , end = "" )
84
94
for nlambda in nlambda_list :
85
95
for ndigit in ndigit_list :
86
96
b = bases [(nlambda , ndigit )]
87
97
sig = f"nlambda{ nlambda } _ndigit{ ndigit } "
88
98
print (
89
99
f"""\
90
- double precision :: s_{ sig } ({ b .size } )
91
- double precision :: tau_{ sig } ({ b .ntau } )
92
- integer :: freq_f_{ sig } ({ b .nfreq_f } )
93
- integer :: freq_b_{ sig } ({ b .nfreq_b } )
94
- double precision :: omega_{ sig } ({ b .nomega } )
95
- double precision :: u_r_{ sig } ({ b .ntau_reduced } * { b .size } )
96
- double precision :: uhat_f_r_{ sig } ({ b .nfreq_f_reduced } * { b .size } )
97
- double precision :: uhat_b_r_{ sig } ({ b .nfreq_b_reduced } * { b .size } )
98
- double precision :: v_r_{ sig } ({ b .nomega_reduced } * { b .size } )
99
- """
100
- )
100
+ REAL(KIND = DP) :: s_{ sig } ({ b .size } )
101
+ REAL(KIND = DP) :: tau_{ sig } ({ b .ntau } )
102
+ INTEGER :: freq_f_{ sig } ({ b .nfreq_f } )
103
+ INTEGER :: freq_b_{ sig } ({ b .nfreq_b } )
104
+ REAL(KIND = DP) :: omega_{ sig } ({ b .nomega } )
105
+ REAL(KIND = DP) :: u_r_{ sig } ({ b .ntau_reduced } * { b .size } )
106
+ REAL(KIND = DP) :: uhat_f_r_{ sig } ({ b .nfreq_f_reduced } * { b .size } )
107
+ REAL(KIND = DP) :: uhat_b_r_{ sig } ({ b .nfreq_b_reduced } * { b .size } )
108
+ REAL(KIND = DP) :: v_r_{ sig } ({ b .nomega_reduced } * { b .size } )
109
+ """ , end = "" )
101
110
102
111
print (
103
112
"""\
104
- contains
105
-
106
- function mk_ir_preset(nlambda, ndigit, beta) result(obj)
107
- integer, intent(in) :: nlambda, ndigit
108
- double precision, intent(in) :: beta
109
- type(IR) :: obj
110
- """ )
113
+ !
114
+ CONTAINS
115
+ !
116
+ !-----------------------------------------------------------------------
117
+ FUNCTION mk_ir_preset(nlambda, ndigit, beta, positive_only) result(obj)
118
+ !-----------------------------------------------------------------------
119
+ !
120
+ INTEGER, INTENT(IN) :: nlambda, ndigit
121
+ REAL(KIND = DP), INTENT(IN) :: beta
122
+ LOGICAL, INTENT(IN), OPTIONAL :: positive_only
123
+ TYPE(IR) :: obj
124
+ !
125
+ """ , end = "" )
111
126
112
127
for nlambda in nlambda_list :
113
128
for ndigit in ndigit_list :
114
129
print (
115
130
f"""\
116
- if (nlambda == { nlambda } .and. ndigit == { ndigit } ) then
117
- obj = mk_nlambda{ nlambda } _ndigit{ ndigit } (beta)
118
- return
119
- end if
120
- """
121
- )
131
+ IF (nlambda == { nlambda } .and. ndigit == { ndigit } ) THEN
132
+ IF ((.NOT. PRESENT(positive_only))) THEN
133
+ obj = mk_nlambda{ nlambda } _ndigit{ ndigit } (beta)
134
+ ELSE
135
+ obj = mk_nlambda{ nlambda } _ndigit{ ndigit } (beta, positive_only)
136
+ ENDIF
137
+ RETURN
138
+ ENDIF
139
+ !
140
+ """ , end = "" )
122
141
123
142
124
143
print (
125
144
"""\
126
- stop "Invalid parameters"
127
- end function
128
- """
129
- )
145
+ STOP "Invalid parameters"
146
+ !
147
+ !-----------------------------------------------------------------------
148
+ END FUNCTION mk_ir_preset
149
+ !-----------------------------------------------------------------------
150
+ !
151
+ """ , end = "" )
130
152
131
153
for nlambda in nlambda_list :
132
154
for ndigit in ndigit_list :
133
155
print_data (nlambda , ndigit , bases [(nlambda , ndigit )])
134
156
135
- print ("end" )
157
+ print (
158
+ """\
159
+ !-----------------------------------------------------------------------
160
+ END MODULE sparse_ir_preset
161
+ !-----------------------------------------------------------------------
162
+ """ , end = "" )
136
163
137
164
138
165
def print_data (nlambda , ndigit , b ):
139
166
sig = f"nlambda{ nlambda } _ndigit{ ndigit } "
140
167
print (
141
- f"""
142
- function mk_nlambda{ nlambda } _ndigit{ ndigit } (beta) result(obj)
143
- double precision, intent(in) :: beta
144
- type(IR) :: obj
145
- complex(kind(0d0)), allocatable :: u(:, :), uhat_f(:, :), uhat_b(:, :), v(:, :), dlr(:, :)
146
- integer, parameter :: size = { b .size } , ntau = { b .ntau } , nfreq_f = { b .nfreq_f } , nfreq_b = { b .nfreq_b } , nomega = { b .nomega }
147
- integer, parameter :: nlambda = { nlambda } , ndigit = { ndigit }
148
- integer, parameter :: ntau_reduced = ntau/2+1, nfreq_f_reduced = nfreq_f/2+1, nfreq_b_reduced = nfreq_b/2+1
149
- integer, parameter :: nomega_reduced = nomega/2+1
150
- double precision, parameter :: lambda = 1.d{ nlambda } , eps = 1.d-{ ndigit }
151
-
152
- integer :: itau, l, ifreq, iomega
153
- """ )
168
+ f"""\
169
+ !-----------------------------------------------------------------------
170
+ FUNCTION mk_nlambda{ nlambda } _ndigit{ ndigit } (beta, positive_only) result(obj)
171
+ !-----------------------------------------------------------------------
172
+ !
173
+ REAL(KIND = DP), INTENT(IN) :: beta
174
+ LOGICAL, INTENT(IN), OPTIONAL :: positive_only
175
+ TYPE(IR) :: obj
176
+ REAL(KIND = DP), ALLOCATABLE :: u(:, :), v(:, :), dlr(:, :)
177
+ COMPLEX(KIND = DP), ALLOCATABLE :: uhat_f(:, :), uhat_b(:, :)
178
+ INTEGER, PARAMETER :: size = { b .size } , ntau = { b .ntau } , nfreq_f = { b .nfreq_f } , nfreq_b = { b .nfreq_b } , nomega = { b .nomega }
179
+ INTEGER, PARAMETER :: nlambda = { nlambda } , ndigit = { ndigit }
180
+ INTEGER, PARAMETER :: ntau_reduced = ntau/2+1
181
+ INTEGER, PARAMETER :: nfreq_f_reduced = nfreq_f/2+1, nfreq_b_reduced = nfreq_b/2+1
182
+ INTEGER, PARAMETER :: nomega_reduced = nomega/2+1
183
+ REAL(KIND = DP), PARAMETER :: lambda = 1.d{ nlambda } , eps = 1.d-{ ndigit }
184
+ !
185
+ INTEGER :: itau, l, ifreq, iomega
186
+ !
187
+ """ , end = "" )
154
188
for varname in ["s" , "tau" , "freq_f" , "freq_b" , "u_r" , "uhat_f_r" , "uhat_b_r" , "omega" , "v_r" ]:
155
- print (8 * " " + f"call init_{ varname } _{ sig } ()" )
189
+ print (2 * " " + f"CALL init_{ varname } _{ sig } ()" )
156
190
157
191
print (
158
- f"""
159
- allocate(u(ntau, size))
160
- allocate(uhat_f(nfreq_f, size))
161
- allocate(uhat_b(nfreq_b, size))
162
- allocate(v(nomega, size))
163
- allocate(dlr(nomega, size))
164
-
165
- ! Use the fact U_l(tau) is even/odd for even/odd l-1.
166
- do l = 1, size
167
- do itau = 1, ntau_reduced
168
- u(itau, l) = u_r_{ sig } (itau + { b .ntau_reduced } *(l-1))
169
- u(ntau-itau+1, l) = (-1)**(l-1) * u_r_{ sig } (itau + { b .ntau_reduced } *(l-1))
170
- end do
171
- end do
172
-
173
- ! Use the fact U^F_l(iv) is pure imaginary/real for even/odd l-1.
174
- do l = 1, size, 2
175
- do ifreq = 1, nfreq_f_reduced
176
- uhat_f(ifreq, l) = cmplx(0.0, uhat_f_r_{ sig } (ifreq + { b .nfreq_f_reduced } *(l-1)), kind(0d0))
177
- end do
178
- end do
179
- do l = 2, size, 2
180
- do ifreq = 1, nfreq_f_reduced
181
- uhat_f(ifreq, l) = cmplx(uhat_f_r_{ sig } (ifreq + { b .nfreq_f_reduced } *(l-1)), 0.0, kind(0d0))
182
- end do
183
- end do
184
- do l = 1, size
185
- do ifreq = 1, nfreq_f
186
- uhat_f(nfreq_f-ifreq+1, l) = conjg(uhat_f(ifreq, l))
187
- end do
188
- end do
189
-
190
- ! Use the fact U^B_l(iv) is pure real/imaginary for even/odd l-1
191
- do l = 1, size, 2
192
- do ifreq = 1, nfreq_b_reduced
193
- uhat_b(ifreq, l) = cmplx(uhat_b_r_{ sig } (ifreq + { b .nfreq_b_reduced } *(l-1)), 0.0d0, kind(0d0))
194
- end do
195
- end do
196
- do l = 2, size, 2
197
- do ifreq = 1, nfreq_b_reduced
198
- uhat_b(ifreq, l) = cmplx(0.0d0, uhat_b_r_{ sig } (ifreq + { b .nfreq_b_reduced } *(l-1)), kind(0d0))
199
- end do
200
- end do
201
- do l = 1, size
202
- do ifreq = 1, nfreq_b
203
- uhat_b(nfreq_b-ifreq+1, l) = conjg(uhat_b(ifreq, l))
204
- end do
205
- end do
206
-
207
- ! Use the fact V_l(omega) is even/odd for even/odd l-1.
208
- do l = 1, size
209
- do iomega = 1, nomega_reduced
210
- v(iomega, l) = v_r_{ sig } (iomega + { b .nomega_reduced } *(l-1))
211
- v(nomega-iomega+1, l) = (-1)**(l-1) * v_r_{ sig } (iomega + { b .nomega_reduced } *(l-1))
212
- end do
213
- do iomega = 1, nomega
214
- dlr(iomega, l) = - s_{ sig } (l) * v(iomega, l)
215
- end do
216
- end do
217
-
218
- call init_ir(obj, beta, lambda, eps,&
219
- s_{ sig } , tau_{ sig } ,&
220
- freq_f_{ sig } , freq_b_{ sig } ,&
221
- u, uhat_f, uhat_b, omega_{ sig } ,&
222
- v, dlr, 1d-20)
223
-
224
- deallocate(u, uhat_f, uhat_b, v, dlr)
225
- end function
226
- """
227
- )
192
+ f"""\
193
+ ALLOCATE(u(ntau, size))
194
+ ALLOCATE(uhat_f(nfreq_f, size))
195
+ ALLOCATE(uhat_b(nfreq_b, size))
196
+ ALLOCATE(v(nomega, size))
197
+ ALLOCATE(dlr(nomega, size))
198
+ !
199
+ ! Use the fact U_l(tau) is even/odd for even/odd l-1.
200
+ DO l = 1, size
201
+ DO itau = 1, ntau_reduced
202
+ u(itau, l) = u_r_{ sig } (itau + { b .ntau_reduced } *(l-1))
203
+ u(ntau-itau+1, l) = (-1)**(l-1) * u_r_{ sig } (itau + { b .ntau_reduced } *(l-1))
204
+ ENDDO
205
+ ENDDO
206
+ !
207
+ ! Use the fact U^F_l(iv) is pure imaginary/real for even/odd l-1.
208
+ DO l = 1, size, 2
209
+ DO ifreq = 1, nfreq_f_reduced
210
+ uhat_f(ifreq, l) = CMPLX(0.0d0, uhat_f_r_{ sig } (ifreq + { b .nfreq_f_reduced } *(l-1)), KIND = DP)
211
+ ENDDO
212
+ ENDDO
213
+ DO l = 2, size, 2
214
+ DO ifreq = 1, nfreq_f_reduced
215
+ uhat_f(ifreq, l) = CMPLX(uhat_f_r_{ sig } (ifreq + { b .nfreq_f_reduced } *(l-1)), 0.0, KIND = DP)
216
+ ENDDO
217
+ ENDDO
218
+ DO l = 1, size
219
+ DO ifreq = 1, nfreq_f
220
+ uhat_f(nfreq_f-ifreq+1, l) = CONJG(uhat_f(ifreq, l))
221
+ ENDDO
222
+ ENDDO
223
+ !
224
+ ! Use the fact U^B_l(iv) is pure real/imaginary for even/odd l-1
225
+ DO l = 1, size, 2
226
+ DO ifreq = 1, nfreq_b_reduced
227
+ uhat_b(ifreq, l) = CMPLX(uhat_b_r_{ sig } (ifreq + { b .nfreq_b_reduced } *(l-1)), 0.0d0, KIND = DP)
228
+ ENDDO
229
+ ENDDO
230
+ DO l = 2, size, 2
231
+ DO ifreq = 1, nfreq_b_reduced
232
+ uhat_b(ifreq, l) = CMPLX(0.0d0, uhat_b_r_{ sig } (ifreq + { b .nfreq_b_reduced } *(l-1)), KIND = DP)
233
+ ENDDO
234
+ ENDDO
235
+ DO l = 1, size
236
+ DO ifreq = 1, nfreq_b
237
+ uhat_b(nfreq_b-ifreq+1, l) = CONJG(uhat_b(ifreq, l))
238
+ ENDDO
239
+ ENDDO
240
+ !
241
+ ! Use the fact V_l(omega) is even/odd for even/odd l-1.
242
+ DO l = 1, size
243
+ DO iomega = 1, nomega_reduced
244
+ v(iomega, l) = v_r_{ sig } (iomega + { b .nomega_reduced } *(l-1))
245
+ v(nomega-iomega+1, l) = (-1)**(l-1) * v_r_{ sig } (iomega + { b .nomega_reduced } *(l-1))
246
+ ENDDO
247
+ DO iomega = 1, nomega
248
+ dlr(iomega, l) = - s_{ sig } (l) * v(iomega, l)
249
+ ENDDO
250
+ ENDDO
251
+ !
252
+ IF ((.NOT. PRESENT(positive_only))) THEN
253
+ CALL init_ir(obj, beta, lambda, eps,&
254
+ s_{ sig } , tau_{ sig } ,&
255
+ freq_f_{ sig } , freq_b_{ sig } ,&
256
+ u, uhat_f, uhat_b, omega_{ sig } ,&
257
+ v, dlr, 1d-20)
258
+ ELSE
259
+ CALL init_ir(obj, beta, lambda, eps,&
260
+ s_{ sig } , tau_{ sig } ,&
261
+ freq_f_{ sig } , freq_b_{ sig } ,&
262
+ u, uhat_f, uhat_b, omega_{ sig } ,&
263
+ v, dlr, 1d-20, positive_only)
264
+ ENDIF
265
+ !
266
+ DEALLOCATE(u, uhat_f, uhat_b, v, dlr)
267
+ !
268
+ !-----------------------------------------------------------------------
269
+ END FUNCTION mk_nlambda{ nlambda } _ndigit{ ndigit }
270
+ !-----------------------------------------------------------------------
271
+ !
272
+ """ , end = "" )
228
273
229
274
print_vector_data (b .s , f"s_{ sig } " )
230
275
print_vector_data (b .x , f"tau_{ sig } " )
@@ -253,25 +298,30 @@ def print_vector_data(vec, var_name):
253
298
vec = np .asfortranarray (vec )
254
299
vec = vec .T .ravel ()
255
300
print (
256
- f"""
257
- subroutine init_{ var_name } ()
258
- """
259
- )
301
+ f"""\
302
+ !-----------------------------------------------------------------------
303
+ SUBROUTINE init_{ var_name } ()
304
+ !-----------------------------------------------------------------------
305
+ !
306
+ """ , end = "" )
260
307
start = 0
261
308
while start < vec .size :
262
309
end = min (start + NLINE , vec .size )
263
310
sub_vec = vec [start :end ]
264
- print (2 * 4 * " " , f"{ var_name } ({ start + 1 } :{ end } ) = (/ &" )
311
+ print (2 * " " + f"{ var_name } ({ start + 1 } :{ end } ) = (/ &" )
265
312
end_str = "&"
266
313
print (", &\n " .join (_array_to_strings (sub_vec )) + end_str )
267
- print (2 * 4 * " " + "/)" )
314
+ print (2 * " " + "/)" )
268
315
start = end
269
316
270
317
print (
271
- f"""
272
- end subroutine
273
- """
274
- )
318
+ f"""\
319
+ !
320
+ !-----------------------------------------------------------------------
321
+ END SUBROUTINE init_{ var_name }
322
+ !-----------------------------------------------------------------------
323
+ !
324
+ """ , end = "" )
275
325
276
326
277
327
if __name__ == '__main__' :
0 commit comments