Choral
ode_solution_mod.f90
Go to the documentation of this file.
1 
56 
58 
59  use real_type
60  use basic_tools
61  use ode_def
62  use ode_problem_mod
63 
64  implicit none
65  private
66 
67  public :: ode_solution
68  public :: clear, print
69 
70  public :: ode_lin_solver
71  public :: ode_nl_1s_solver
72  public :: ode_nl_ms_solver
75 
76  ! %----------------------------------------%
77  ! | |
78  ! | DERIVED TYPE |
79  ! | |
80  ! %----------------------------------------%
81 
86  type :: ode_solution
87 
89  integer :: type = -1
90 
92  integer :: n = -1
93 
95  integer :: na = -1
96 
99  integer :: dof = -1
100 
102  real(RP), dimension(:,:,:), allocatable :: y
103 
105  real(RP), dimension(:,:,:), allocatable :: ay
106 
108  real(RP), dimension(:,:,:), allocatable :: by
109 
111  real(RP), dimension(:,:) , allocatable :: v
112 
114  real(RP), dimension(:), allocatable :: rhs
115 
117  real(RP), dimension(:), allocatable :: aux
118 
120  integer :: nv=0, ny=0, nfy=0
121 
127  integer :: ierr = 0
128 
132  integer, dimension(:), allocatable :: y_i
133 
137  integer, dimension(:), allocatable :: f_i
138 
142  integer, dimension(:), allocatable :: v_i
143 
144  contains
145 
148 
149  end type ode_solution
150 
151 
152  ! %----------------------------------------%
153  ! | |
154  ! | GENERIc SUBROUTINES |
155  ! | |
156  ! %----------------------------------------%
157  interface clear
158  module procedure ode_solution_clear
159  end interface clear
160 
161  interface ode_solution
162  module procedure ode_solution_create
163  end interface ode_solution
164 
165  interface print
166  module procedure ode_solution_print
167  end interface print
168 
169 
170  ! %----------------------------------------%
171  ! | |
172  ! | ABSTRACT INTERFACES |
173  ! | |
174  ! %----------------------------------------%
175  abstract interface
176 
177  subroutine ode_output_proc(tn, s, stop)
178  import :: rp, ode_solution
179  real(RP) , intent(in) :: tn
180  type(ode_solution), intent(in) :: s
181  logical , intent(inout) :: stop
182  end subroutine ode_output_proc
183 
184  subroutine ode_lin_solver(sol, ierr, dt, pb, KInv)
185  import :: rp, ode_problem, ode_solution, linsystem_solver
186  type(ode_solution) , intent(inout) :: sol
187  logical , intent(out) :: ierr
188  real(RP) , intent(in) :: dt
189  type(ode_problem) , intent(in) :: pb
190  procedure(linsystem_solver) :: KInv
191  end subroutine ode_lin_solver
192 
193  subroutine ode_nl_1s_solver(sol, dt, t, pb)
195  type(ode_solution), intent(inout) :: sol
196  real(RP) , intent(in) :: dt, t
197  type(ode_problem) , intent(in) :: pb
198  end subroutine ode_nl_1s_solver
199 
201  subroutine ode_nl_ms_solver(sol, dt, N, Na)
202  import :: rp, ode_solution
203  type(ode_solution), intent(inout) :: sol
204  real(RP) , intent(in) :: dt
205  integer , intent(in) :: N, Na
206  end subroutine ode_nl_ms_solver
207 
208  end interface
209 
210 contains
211 
214  subroutine void_ode_output(tn, s, stop)
215  real(RP) , intent(in) :: tn
216  type(ode_solution), intent(in) :: s
217  logical , intent(inout) :: stop
218 
219  stop = .false.
220  end subroutine void_ode_output
221 
222 
225  subroutine ode_solution_clear(sol)
226  type(ode_solution), intent(inout) :: sol
227 
228  call freemem( sol%Y )
229  call freemem( sol%AY )
230  call freemem( sol%BY )
231  call freemem( sol%V )
232  call freemem( sol%aux )
233  call freemem( sol%rhs )
234 
235  call freemem( sol%Y_i )
236  call freemem( sol%F_i )
237  call freemem( sol%V_i )
238 
239  sol%N = -1
240  sol%Na = -1
241  sol%dof = -1
242 
243  sol%nV = 0
244  sol%nY = 0
245  sol%nFY = 0
246 
247  sol%ierr = 0
248 
249  end subroutine ode_solution_clear
250 
251 
257  function ode_solution_create(pb, nV, NY, NFY) result(sol)
258  type(ode_solution) :: sol
259  type(ode_problem) , intent(in) :: pb
260  integer, intent(in), optional :: nV, NY, NFY
261 
262  call clear(sol)
263 
264  sol%type = pb%type
265  sol%N = pb%N
266  sol%Na = pb%Na
267  sol%dof = pb%dof
268 
269  if (present(nv)) then
270  if ( nv < 0) call quit(&
271  & "ode_solution_mod: ode_solution_create:&
272  & parameter 'nV' must be >= 0" )
273 
274  sol%nV = nv
275 
276  call allocmem(sol%V , pb%dof, nv)
277  call allocmem(sol%aux, pb%dof)
278  call allocmem(sol%rhs, pb%dof)
279  call allocmem(sol%V_i, nv)
280 
281  end if
282 
283  if (present(ny)) then
284 
285  if ( ny < 0) call quit(&
286  & "ode_solution_mod: ode_solution_create::&
287  & parameter 'nY' must be >= 0" )
288 
289  sol%nY = ny
290  call allocmem(sol%Y, pb%N, ny, pb%dof)
291  call allocmem(sol%Y_i, ny)
292 
293  end if
294 
295  if (present(nfy)) then
296 
297  if ( nfy < 0) call quit(&
298  & "ode_solution_mod: ode_solution_create::&
299  & parameter 'nFY' must be >= 0" )
300 
301  sol%nFY = nfy
302  call allocmem(sol%BY, pb%N , nfy, pb%dof)
303  call allocmem(sol%AY, pb%Na, nfy, pb%dof)
304  call allocmem(sol%F_i, nfy)
305 
306  end if
307 
308  call ode_solution_init_indexes(sol)
309 
310  end function ode_solution_create
311 
312 
315  subroutine ode_solution_init_indexes(sol)
316  type(ode_solution) , intent(inout) :: sol
317 
318  integer :: ii
319 
320  if (sol%nY > 0) then
321  do ii=1, sol%nY
322  sol%Y_i(ii) = ii
323  end do
324  end if
325 
326  if (sol%nFY > 0) then
327  do ii=1, sol%nFY
328  sol%F_i(ii) = ii
329  end do
330  end if
331 
332  if (sol%nV > 0) then
333  do ii=1, sol%nV
334  sol%V_i(ii) = ii
335  end do
336  end if
337 
338 
339  end subroutine ode_solution_init_indexes
340 
343  subroutine ode_solution_print(sol)
344  type(ode_solution) , intent(in) :: sol
345 
346  write(*,*)"ode_solution_mod: ode_solution_print"
347 
348  if (allocated(sol%Y)) then
349  write(*,*)" Array Y , shape =", &
350  & shape(sol%Y)
351  end if
352 
353  if (allocated(sol%AY)) then
354  write(*,*)" Array AY, shape =", &
355  & shape(sol%AY)
356  end if
357 
358  if (allocated(sol%BY)) then
359  write(*,*)" Array BY, shape =", &
360  & shape(sol%BY)
361  end if
362 
363 
364  if (allocated(sol%V)) then
365  write(*,*)" Array V , shape =", &
366  & " ", shape(sol%V)
367  end if
368 
369  end subroutine ode_solution_print
370 
371 end module ode_solution_mod
Type ode_solution: data structure to solve ODE/PDE problems.
BASIC TOOLS
Definition: basic_tools.F90:8
deallocate memory for real(RP) arrays
Definition: real_type.F90:163
DERIVED TYPE ode_problem: definition of ODE/PDE problems
integer, parameter, public rp
real(kind=RP) = real precision in the code REAL_TOL = epsilon to test real equality ...
Definition: real_type.F90:90
subroutine, public quit(message)
Stop program execution, display an error messahe.
REAL NUMBERS PRECISION IN CHORAL: selects simple/double/quad
Definition: real_type.F90:33
subroutine, public ode_solution_init_indexes(sol)
initialise the ode_indexes
DERIVED TYPE ode_solution: data straucture to solve ODEs
type(ode_solution) function ode_solution_create(pb, nV, NY, NFY)
Bottom level constructor for ode_solution.
subroutine ode_solution_clear(sol)
destructor
allocate memory for real(RP) arrays
Definition: real_type.F90:158
BOTTOM LEVEL MODULE FOR ODEs
Definition: ode_def.f90:12
subroutine ode_solution_print(sol)
print ode_solution
Type ode_problem: definition of ODE/PDE problems.
subroutine, public void_ode_output(tn, s, stop)
void output for ode resolution