This repository has been archived by the owner on Jul 20, 2021. It is now read-only.
-
Notifications
You must be signed in to change notification settings - Fork 5
/
Copy pathvolume_averages.f90
121 lines (83 loc) · 2.85 KB
/
volume_averages.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
module volume_averages
public :: init_volume_averages, finish_volume_averages
public :: fieldline_average
public :: volume_average
! public :: flux_surface_average
public :: mode_fac
private
interface fieldline_average
module procedure fieldline_average_real
module procedure fieldline_average_complex
end interface
real, dimension (:), allocatable :: mode_fac
contains
subroutine init_volume_averages
use kt_grids, only: aky, naky
implicit none
if (.not.allocated(mode_fac)) then
allocate (mode_fac(naky)) ; mode_fac = 2.0
if (aky(1)<epsilon(0.)) mode_fac(1) = 1.0
end if
end subroutine init_volume_averages
subroutine finish_volume_averages
implicit none
if (allocated(mode_fac)) deallocate (mode_fac)
end subroutine finish_volume_averages
!==============================================
!============ FIELD LINE AVERAGE ==============
!==============================================
subroutine fieldline_average_real (unavg, avg)
use zgrid, only: nzgrid, ntubes
use kt_grids, only: nakx, naky
use stella_geometry, only: dl_over_b
implicit none
real, dimension (:,:,-nzgrid:,:), intent (in) :: unavg
real, dimension (:,:), intent (out) :: avg
integer :: it, ia
ia = 1
avg = 0.0
do it = 1, ntubes
avg = avg + sum(spread(spread(dl_over_b(ia,:),1,naky),2,nakx)*unavg(:,:,:,it),dim=3)
end do
avg = avg/real(ntubes)
end subroutine fieldline_average_real
subroutine fieldline_average_complex (unavg, avg)
use zgrid, only: nzgrid, ntubes
use kt_grids, only: nakx, naky
use stella_geometry, only: dl_over_b
implicit none
complex, dimension (:,:,-nzgrid:,:), intent (in) :: unavg
complex, dimension (:,:), intent (out) :: avg
integer :: it, ia
ia = 1
avg = 0.0
do it = 1, ntubes
avg = avg + sum(spread(spread(dl_over_b(ia,:),1,naky),2,nakx)*unavg(:,:,:,it),dim=3)
end do
avg = avg/real(ntubes)
end subroutine fieldline_average_complex
!==============================================
!============== VOLUME AVERAGE ================
!==============================================
subroutine volume_average (unavg, avg)
use zgrid, only: nzgrid, ntubes
use kt_grids, only: naky, nakx
use stella_geometry, only: dl_over_b
implicit none
complex, dimension (:,:,-nzgrid:,:), intent (in) :: unavg
real, intent (out) :: avg
integer :: iky, ikx, iz, it, ia
ia = 1
avg = 0.
do it = 1, ntubes
do iz = -nzgrid, nzgrid
do ikx = 1, nakx
do iky = 1, naky
avg = avg + real(unavg(iky,ikx,iz,it)*conjg(unavg(iky,ikx,iz,it)))*mode_fac(iky)*dl_over_b(ia,iz)
end do
end do
end do
end do
avg = avg/real(ntubes)
end subroutine volume_average
end module volume_averages