Report Typos and Errors    
Semi-Lagrangian Library
Modular library for kinetic and gyrokinetic simulations of plasmas in fusion energy devices.
sll_m_preconditioner_fft.F90
Go to the documentation of this file.
1 
9 
11  !+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
12 #include "sll_assert.h"
13 #include "sll_memory.h"
14 #include "sll_working_precision.h"
15 
16 
19 
20  use sll_m_spline_fem_utilities, only : &
23 
24  use sll_m_linear_solver_block, only : &
26 
27  implicit none
28 
29  public :: &
31 
32  private
33  !+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
34 
36  sll_real64 :: delta_x(3)
37  sll_int32 :: n_dofs0(3)
38  sll_int32 :: n_dofs1(3)
39  sll_int32 :: s_deg_0(3)
40  sll_int32 :: s_deg_1(3)
41  type(sll_t_linear_solver_spline_mass_fft) :: inverse_mass1_1d(3)
42  type(sll_t_linear_solver_spline_mass_fft) :: inverse_mass2_1d(3)
43  type(sll_t_linear_solver_block) :: inverse_mass1_3d
44  type(sll_t_linear_solver_block) :: inverse_mass2_3d
45  sll_real64, allocatable :: mass_line0_1(:)
46  sll_real64, allocatable :: mass_line0_2(:)
47  sll_real64, allocatable :: mass_line0_3(:)
48  sll_real64, allocatable :: mass_line1_1(:)
49  sll_real64, allocatable :: mass_line1_2(:)
50  sll_real64, allocatable :: mass_line1_3(:)
51 
52  contains
53 
54  procedure :: &
55  init => init_3d_trafo
56  procedure :: &
57  free => free_3d_trafo
58 
59 
61 
62 contains
63 
64 
65  subroutine init_3d_trafo( self, Lx, n_cells, s_deg_0, clamped )
66  class(sll_t_preconditioner_fft), intent( inout ) :: self
67  sll_real64, intent( in ) :: lx(3)
68  sll_int32, intent( in ) :: n_cells(3)
69  sll_int32, intent( in ) :: s_deg_0(3)
70  logical, optional :: clamped
71  ! local variables
72  sll_int32 :: j
73 
74  sll_real64, allocatable :: eig_values_mass_0_1(:)
75  sll_real64, allocatable :: eig_values_mass_0_2(:)
76  sll_real64, allocatable :: eig_values_mass_0_3(:)
77  sll_real64, allocatable :: eig_values_mass_1_1(:)
78  sll_real64, allocatable :: eig_values_mass_1_2(:)
79  sll_real64, allocatable :: eig_values_mass_1_3(:)
80 
81 
82  self%delta_x = lx / real(n_cells,f64)
83  self%s_deg_0 = s_deg_0
84  self%s_deg_1 = s_deg_0 - 1
85  self%n_dofs0 = n_cells
86  self%n_dofs1 = n_cells
87 
88  if( present(clamped) )then
89  if(clamped) then
90  self%n_dofs0(1) = n_cells(1)+self%s_deg_0(1)
91  self%n_dofs1(1) = n_cells(1)+self%s_deg_1(1)
92  end if
93  end if
94  ! Allocate scratch data
95  allocate( self%mass_line0_1(s_deg_0(1)+1) )
96  allocate( self%mass_line1_1(s_deg_0(1)) )
97  allocate( self%mass_line0_2(s_deg_0(2)+1) )
98  allocate( self%mass_line1_2(s_deg_0(2)) )
99  allocate( self%mass_line0_3(s_deg_0(3)+1) )
100  allocate( self%mass_line1_3(s_deg_0(3)) )
101  call sll_s_spline_fem_mass_line ( self%s_deg_0(1), self%mass_line0_1 )
102  call sll_s_spline_fem_mass_line ( self%s_deg_1(1), self%mass_line1_1 )
103 
104  call sll_s_spline_fem_mass_line ( self%s_deg_0(2), self%mass_line0_2 )
105  call sll_s_spline_fem_mass_line ( self%s_deg_1(2), self%mass_line1_2 )
106 
107  call sll_s_spline_fem_mass_line ( self%s_deg_0(3), self%mass_line0_3 )
108  call sll_s_spline_fem_mass_line ( self%s_deg_1(3), self%mass_line1_3 )
109 
110  self%mass_line0_1 = self%delta_x(1) * self%mass_line0_1
111  self%mass_line1_1 = self%delta_x(1) * self%mass_line1_1
112 
113  self%mass_line0_2 = self%delta_x(2) * self%mass_line0_2
114  self%mass_line1_2 = self%delta_x(2) * self%mass_line1_2
115 
116  self%mass_line0_3 = self%delta_x(3) * self%mass_line0_3
117  self%mass_line1_3 = self%delta_x(3) * self%mass_line1_3
118 
119  allocate( eig_values_mass_0_1( self%n_dofs0(1) ) )
120  allocate( eig_values_mass_0_2( self%n_dofs0(2) ) )
121  allocate( eig_values_mass_0_3( self%n_dofs0(3) ) )
122  allocate( eig_values_mass_1_1( self%n_dofs1(1) ) )
123  allocate( eig_values_mass_1_2( self%n_dofs1(2) ) )
124  allocate( eig_values_mass_1_3( self%n_dofs1(3) ) )
125 
126  call sll_s_spline_fem_compute_mass_eig( self%n_dofs0(1), self%s_deg_0(1), self%mass_line0_1, &
127  eig_values_mass_0_1 )
128  call sll_s_spline_fem_compute_mass_eig( self%n_dofs0(2), self%s_deg_0(2), self%mass_line0_2, &
129  eig_values_mass_0_2 )
130  call sll_s_spline_fem_compute_mass_eig( self%n_dofs0(3), self%s_deg_0(3), self%mass_line0_3, &
131  eig_values_mass_0_3 )
132  call sll_s_spline_fem_compute_mass_eig( self%n_dofs1(1), self%s_deg_1(1), self%mass_line1_1, &
133  eig_values_mass_1_1 )
134  call sll_s_spline_fem_compute_mass_eig( self%n_dofs1(2), self%s_deg_1(2), self%mass_line1_2, &
135  eig_values_mass_1_2 )
136  call sll_s_spline_fem_compute_mass_eig( self%n_dofs1(3), self%s_deg_1(3), self%mass_line1_3, &
137  eig_values_mass_1_3 )
138 
139  eig_values_mass_0_1 = 1._f64/eig_values_mass_0_1
140  eig_values_mass_0_2 = 1._f64/eig_values_mass_0_2
141  eig_values_mass_0_3 = 1._f64/eig_values_mass_0_3
142  eig_values_mass_1_1 = 1._f64/eig_values_mass_1_1
143  eig_values_mass_1_2 = 1._f64/eig_values_mass_1_2
144  eig_values_mass_1_3 = 1._f64/eig_values_mass_1_3
145 
146  call self%inverse_mass1_1d(1)%create( self%n_dofs1, eig_values_mass_1_1, eig_values_mass_0_2, eig_values_mass_0_3 )
147  call self%inverse_mass1_1d(2)%create( self%n_dofs0, eig_values_mass_0_1, eig_values_mass_1_2, eig_values_mass_0_3 )
148  call self%inverse_mass1_1d(3)%create( self%n_dofs0, eig_values_mass_0_1, eig_values_mass_0_2, eig_values_mass_1_3 )
149 
150  call self%inverse_mass2_1d(1)%create( self%n_dofs0, eig_values_mass_0_1, eig_values_mass_1_2, eig_values_mass_1_3 )
151  call self%inverse_mass2_1d(2)%create( self%n_dofs1, eig_values_mass_1_1, eig_values_mass_0_2, eig_values_mass_1_3 )
152  call self%inverse_mass2_1d(3)%create( self%n_dofs1, eig_values_mass_1_1, eig_values_mass_1_2, eig_values_mass_0_3 )
153 
154  call self%inverse_mass1_3d%create( 3, 3 )
155  call self%inverse_mass2_3d%create( 3, 3 )
156  do j=1,3
157  call self%inverse_mass1_3d%set( j, j, self%inverse_mass1_1d(j) )
158  call self%inverse_mass2_3d%set( j, j, self%inverse_mass2_1d(j) )
159  end do
160 
161 
162  end subroutine init_3d_trafo
163 
164 
165 
166  subroutine free_3d_trafo( self )
167  class(sll_t_preconditioner_fft) :: self
168  !local variable
169  sll_int32 :: j
170 
171  do j=1, 3
172  call self%inverse_mass1_1d(j)%free()
173  call self%inverse_mass2_1d(j)%free()
174  end do
175 
176  call self%inverse_mass1_3d%free()
177  call self%inverse_mass2_3d%free()
178 
179  end subroutine free_3d_trafo
180 
181 
182 end module sll_m_preconditioner_fft
module for a block linear solver
Invert a circulant matrix based on diagonalization in Fourier space (3d version)
Module interface to solve Maxwell's equations with coordinate transformation in 3D The linear systems...
subroutine init_3d_trafo(self, Lx, n_cells, s_deg_0, clamped)
Utilites for Maxwell solver's with spline finite elements.
subroutine, public sll_s_spline_fem_mass_line(degree, mass_line)
Computes the mass line for a mass matrix with degree splines.
subroutine, public sll_s_spline_fem_compute_mass_eig(n_cells, degree, mass_line, eig_values_mass)
Compute eigenvalues of mass matrix with given line mass_line.
Linear solver for FFT-based inversion of 3d tensor product of circulant matrices (extending the abstr...
    Report Typos and Errors