-
Notifications
You must be signed in to change notification settings - Fork 1
/
Copy pathwpush1.f90
244 lines (244 loc) · 7.82 KB
/
wpush1.f90
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
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
!-----------------------------------------------------------------------
!
module wpush1
!
! Fortran90 wrappers to 1d PIC library push1.f
! wdistr1 calculates initial particle co-ordinates with uniform plasma
! maxwellian velocities
! calls DISTR1
! wgpush1 push particles with leap-frog
! calls GPUSH1L
! wgmpush1 push particles with leap-frog using OpenMP
! calls MGPUSH1L
! wgpost1 deposit charge
! calls GPOST1L
! wcguard1 copy guard cells for periodic 1d scalar data
! calls CGUARD1L
! waguard1 add guard cells for periodic 1d scalar data
! calls AGUARD1L
! pois1_init calculates table needed by 1d poisson solver
! calls POIS1
! wpois1 poisson solver for periodic 1d electric field
! calls POIS1
! fft1_init calculates tables needed by 1d FFTs
! calls WFFT1RINIT
! wfft1r wrapper function for scalar 1d real/complex FFT
! calls FFT1RXX
! meaddext1 add external traveling wave field to electric field for 1d
! code
! calls EADDEXT1
! fprecision determines if default reals are actually doubles
! written by viktor k. decyk, ucla
! copyright 2021, regents of the university of california
! update: december 18, 2023
!
use push1_h
implicit none
!
! t = scratch array for wfft1r
complex, dimension(:), allocatable :: t
integer :: szt = 0
save
!
private :: t, szt
!
contains
!
!-----------------------------------------------------------------------
subroutine wdistr1(part,vtx,vdx,npx,nx)
! calculates initial particle co-ordinates with uniform plasma and
! maxwellian velocities
implicit none
integer, intent(in) :: npx, nx
real, intent(in) :: vtx, vdx
real, dimension(:,:), intent(inout) :: part
! local data
integer :: ipbc = 1
integer :: idimp, nop
! extract dimensions
idimp = size(part,1); nop = size(part,2)
! call low level procedure
call DISTR1(part,vtx,vdx,npx,idimp,nop,nx,ipbc)
end subroutine
!
!-----------------------------------------------------------------------
subroutine wgpush1(part,fx,qbm,dt,ek,nx)
! push particles with leap-frog
implicit none
integer, intent(in) :: nx
real, intent(in) :: qbm, dt
real, dimension(1), intent(inout) :: ek
real, dimension(:,:), intent(inout) :: part
real, dimension(:), intent(in) :: fx
! local data
integer :: ipbc = 1
integer :: idimp, nop, nxv
! extract dimensions
idimp = size(part,1); nop = size(part,2)
nxv = size(fx,1)
! call low level procedure
call GPUSH1L(part,fx,qbm,dt,ek(1),idimp,nop,nx,nxv,ipbc)
end subroutine
!
!-----------------------------------------------------------------------
subroutine wgmpush1(part,fx,qbm,dt,ek,nvp,nx)
! push particles with leap-frog
implicit none
integer, intent(in) :: nvp, nx
real, intent(in) :: qbm, dt
real, dimension(1), intent(inout) :: ek
real, dimension(:,:), intent(inout) :: part
real, dimension(:), intent(in) :: fx
! local data
integer :: ipbc = 1
integer :: idimp, nop, nxv
! extract dimensions
idimp = size(part,1); nop = size(part,2)
nxv = size(fx,1)
! call low level procedure
call MGPUSH1L(part,fx,qbm,dt,ek(1),idimp,nop,nvp,nx,nxv,ipbc)
end subroutine
!
!-----------------------------------------------------------------------
subroutine wgpost1(part,q,qm)
! deposit charge
implicit none
real, intent(in) :: qm
real, dimension(:,:), intent(in) :: part
real, dimension(:), intent(inout) :: q
! local data
integer :: idimp, nop, nxv
! extract dimensions
idimp = size(part,1); nop = size(part,2)
nxv = size(q,1)
! call low level procedure
call GPOST1L(part,q,qm,nop,idimp,nxv)
end subroutine
!
!-----------------------------------------------------------------------
subroutine wcguard1(fx,nx)
! copy guard cells for periodic 1d scalar data
implicit none
integer, intent(in) :: nx
real, dimension(:), intent(inout) :: fx
! local data
integer :: nxe
! extract dimensions
nxe = size(fx,1)
! call low level procedure
call CGUARD1L(fx,nx,nxe)
end subroutine
!
!-----------------------------------------------------------------------
subroutine waguard1(q,nx)
! add guard cells for periodic 1d scalar data
implicit none
integer, intent(in) :: nx
real, dimension(:), intent(inout) :: q
! local data
integer :: nxe
! extract dimensions
nxe = size(q,1)
! call low level procedure
call AGUARD1L(q,nx,nxe)
end subroutine
!
!-----------------------------------------------------------------------
subroutine pois1_init(ffc,ax,affp,nx)
! calculates table needed by 1d poisson solver
implicit none
integer, intent(in) :: nx
real, intent(in) :: ax, affp
complex, dimension(:), intent(inout) :: ffc
! local data
integer :: isign = 0
real, dimension(1) :: we
real, dimension(1) :: q, fx
! call low level procedure
call POIS1(q,fx,isign,ffc,ax,affp,we(1),nx)
end subroutine
!
!-----------------------------------------------------------------------
subroutine wpois1(q,fx,ffc,we,nx)
! poisson solver for periodic 1d electric field
implicit none
integer, intent(in) :: nx
real, dimension(1), intent(inout) :: we
real, dimension(:), intent(in) :: q
real, dimension(:), intent(inout) :: fx
complex, dimension(:), intent(inout) :: ffc
! local data
integer :: isign = -1
real :: ax, affp
! call low level procedure
call POIS1(q,fx,isign,ffc,ax,affp,we(1),nx)
end subroutine
!
!-----------------------------------------------------------------------
subroutine fft1_init(mixup,sct,indx)
! calculates tables needed by 1d FFTs
implicit none
integer, intent(in) :: indx
integer, dimension(:), intent(inout) :: mixup
complex, dimension(:), intent(inout) :: sct
! local data
integer :: nxhd
! extract dimensions
nxhd = size(mixup,1)
! call low level procedure
call WFFT1RINIT(mixup,sct,indx,nxhd)
end subroutine
!
!-----------------------------------------------------------------------
subroutine wfft1r(f,isign,mixup,sct,indx)
! wrapper function for scalar 1d real/complex FFT
implicit none
integer, intent(in) :: isign, indx
real, dimension(:), intent(inout) :: f
integer, dimension(:), intent(in) :: mixup
complex, dimension(:), intent(in) :: sct
! local data
integer :: nxd, nxhd
! extract dimensions
nxd = size(f,1)
nxhd = size(mixup,1)
! check if required size of buffer t has increased
if (szt < nxhd) then
if (szt /= 0) deallocate(t)
! allocate new buffers
allocate(t(nxhd))
szt = nxhd
endif
! call low level procedure
call FFT1RXX(f,t,isign,mixup,sct,indx,nxd,nxhd)
end subroutine
!
!-----------------------------------------------------------------------
subroutine meaddext1(fxe,amodex,freq,time,trmp,toff,el0,er0,nx)
! add external traveling wave field to electric field for 1-2/2d code
implicit none
integer, intent(in) :: nx
real, intent(in) :: amodex, freq, time, trmp, toff, el0, er0
real, dimension(:), intent(inout) :: fxe
! local data
integer :: nxe
! extract dimensions
nxe = size(fxe,1)
! call low level procedure
call EADDEXT1(fxe,amodex,freq,time,trmp,toff,el0,er0,nx,nxe)
end subroutine
!
!-----------------------------------------------------------------------
integer function fprecision()
! function determines if default reals are actually doubles
implicit none
real :: prec
! ndprec = (0,1) = (no,yes) = (normal,autodouble) precision used
if (digits(prec) > 24) then
fprecision = 1
else
fprecision = 0
endif
end function fprecision
!
end module