Neko  0.9.0
A portable framework for high-order spectral element flow simulations
fluid_plan1.f90
Go to the documentation of this file.
1 
3 module fluid_plan1
4  use fluid_scheme
5  implicit none
6 
7  type, extends(fluid_scheme_t) :: fluid_plan1_t
8  type(space_t) :: yh
9  type(dofmap_t) :: dm_yh
10  type(gs_t) :: gs_yh
11  type(coef_t) :: c_yh
13  contains
14  procedure, pass(this) :: init => fluid_plan1_init
15  procedure, pass(this) :: free => fluid_plan1_free
16  procedure, pass(this) :: step => fluid_plan1_step
17  end type fluid_plan1_t
18 
19 contains
20 
21  subroutine fluid_plan1_init(this, msh, lx, param)
22  class(fluid_plan1_t), target, intent(inout) :: this
23  type(mesh_t), target, intent(inout) :: msh
24  integer, intent(inout) :: lx
25  type(param_t), target, intent(inout) :: param
26  character(len=15), parameter :: scheme = 'plan1 (Pn/Pn-2)'
27  integer :: lx2
28 
29  call this%free()
30 
32  call this%scheme_init(msh, lx, param, kspv_init=.true., scheme=scheme)
33 
35  lx2 = lx - 2
36  if (msh%gdim .eq. 2) then
37  call this%Yh%init(gll, lx2, lx2)
38  else
39  call this%Yh%init(gll, lx2, lx2, lx2)
40  end if
41 
42  call this%dm_Yh%init(msh, this%Yh)
43 
44  call this%p%init(this%dm_Yh)
45 
46  call gs_init(this%gs_Yh, this%dm_Yh)
47 
48  call this%c_Yh%init(this%gs_Yh)
49 
50  call fluid_scheme_solver_factory(this%ksp_prs, this%dm_Yh%size(), &
51  param%ksp_prs, param%abstol_prs)
52  call fluid_scheme_precon_factory(this%pc_prs, this%ksp_prs, &
53  this%c_Yh, this%dm_Yh, this%gs_Yh, this%bclst_prs, param%pc_prs)
54 
55 
56  end subroutine fluid_plan1_init
57 
58  subroutine fluid_plan1_free(this)
59  class(fluid_plan1_t), intent(inout) :: this
60 
61  ! Deallocate velocity and pressure fields
62  call this%scheme_free()
63 
64  ! Deallocate Pressure dofmap and space
65  call space_free(this%Yh)
66 
67  call gs_free(this%gs_Yh)
68 
69  call this%c_Yh%free()
70 
71  end subroutine fluid_plan1_free
72 
73  subroutine fluid_plan1_step(this, t, tstep, ext_bdf, dt_controller)
74  class(fluid_plan1_t), target, intent(inout) :: this
75  real(kind=rp), intent(inout) :: t
76  integer, intent(inout) :: tstep
77  type(time_scheme_controller_t), intent(inout) :: ext_bdf
78  type(time_step_controller_t), intent(in) :: dt_controller
79 
80  if (this%freeze) return
81 
82  end subroutine fluid_plan1_step
83 
84 end module fluid_plan1
Classic NEKTON formulation Compute pressure and velocity using consistent approximation spaces.
Definition: fluid_plan1.f90:3
subroutine fluid_plan1_step(this, t, tstep, ext_bdf, dt_controller)
Definition: fluid_plan1.f90:74
subroutine fluid_plan1_init(this, msh, lx, param)
Definition: fluid_plan1.f90:22
subroutine fluid_plan1_free(this)
Definition: fluid_plan1.f90:59
Fluid formulations.
subroutine fluid_scheme_solver_factory(ksp, n, solver, max_iter, abstol, monitor)
Initialize a linear solver.
subroutine fluid_scheme_precon_factory(pc, ksp, coef, dof, gs, bclst, pctype)
Initialize a Krylov preconditioner.
Base type of all fluid formulations.