Choral
ode_solver_mod.f90
Go to the documentation of this file.
1 
106 
108 
109  use choral_constants
110  use choral_variables
111  use real_type
112  use basic_tools
113  use krylov_mod
114  use ode_def
115  use ode_problem_mod
116  use ode_solution_mod
117  use ode_nl_1s_mod
118  use ode_nl_ms_mod
119  use ode_lin_1s_mod
120  use ode_sl_ms_mod
121  use ode_sl_nl_ms_mod
122  use ode_sl_nl_dc_mod
123  use ode_opsplt_mod
124  use ode_output_mod
125 
126  implicit none
127  private
128 
129  ! %----------------------------------------%
130  ! | |
131  ! | PUBLIC DATA |
132  ! | |
133  ! %----------------------------------------%
134 
135  !! DERIVED TYPE
136  public :: ode_solver
137 
138  !! GENERIC SUBROUTINES
139  public :: clear
140  public :: ode_solution
141  public :: valid
142  public :: print
143  public :: solve
144  public :: initialcond
145 
146  !! NON GENERIC SUBROUTINES
147  public :: name_ode_solver_type
148  public :: check_ode_method
149  public :: set_ode_solver_output
150 
151  ! %----------------------------------------%
152  ! | |
153  ! | DERIVED TYPE |
154  ! | |
155  ! %----------------------------------------%
156 
161  type :: ode_solver
162 
164  integer :: verb = 0
165 
167  integer :: type = -1
168 
170  integer :: pb = -1
171 
173  type(ode_opsplt) :: os
174 
176  integer :: l_meth = -1
177 
179  integer :: sl_meth = -1
180 
182  integer :: dc_meth = -1
183 
185  integer :: nl_meth = -1
186 
188  type(krylov) :: kry
189 
191  logical :: check_overflow = .false.
192 
194  procedure(ode_output_proc), nopass, pointer :: output => void_ode_output
195 
196  contains
197 
200 
201  end type ode_solver
202 
203 
204  ! %----------------------------------------%
205  ! | |
206  ! | GENERIc SUBROUTINES |
207  ! | |
208  ! %----------------------------------------%
209 
210  interface clear
211  module procedure ode_solver_clear
212  end interface clear
213 
214  interface ode_solver
215  module procedure ode_solver_create
216  end interface ode_solver
217 
218  interface valid
219  module procedure ode_solver_valid
220  end interface valid
221 
222  interface print
223  module procedure ode_solver_print
224  end interface print
225 
226  interface ode_solution
227  module procedure ode_solution_create
228  end interface ode_solution
229 
231  interface solve
232  module procedure ode_solver_solve
233  end interface solve
234 
236  interface initialcond
237  module procedure homogeneous_initialcond
238  end interface initialcond
239 
240 contains
241 
242 
245  subroutine ode_solver_clear(slv)
246  type(ode_solver), intent(inout) :: slv
247 
248  slv%verb = 0
249  slv%type = -1
250  slv%pb = -1
251  slv%L_meth = -1
252  slv%SL_meth = -1
253  slv%DC_meth = -1
254  slv%NL_meth = -1
255  slv%check_overflow = .false.
256  slv%output => void_ode_output
257  call clear(slv%os)
258  call clear(slv%kry)
259 
260  end subroutine ode_solver_clear
261 
262 
281  function ode_solver_create(pb, type, &
282  & os_meth, SL_meth, NL_meth, L_meth, &
283  & DC_meth, &
284  & check_overflow, verb) result(slv)
285  type(ode_solver) :: slv
286  type(ode_problem), intent(in) :: pb
287  integer , intent(in) :: type
288  integer, optional, intent(in) :: os_meth, SL_meth, DC_meth
289  integer, optional, intent(in) :: NL_meth, L_meth, verb
290  logical, optional, intent(in) :: check_overflow
291 
292  logical :: bool
293 
294  if (choral_verb>0) write(*,*) &
295  &"ode_solver_mod : ode_solver_create"
296 
297  call clear(slv)
298 
299  !! default krylov
300  !!
301  slv%kry = krylov(type=kry_cg, tol=1e-6_rp, itmax=100, verb=0)
302 
303  !! Asociated ODE problem
304  !!
305  bool = valid(pb)
306  if (.NOT.bool) call quit(&
307  & "ode_solver_mod: ode_solver_create: 'pb' not valid")
308  slv%pb = pb%type
309 
310  !! Type
311  !!
312  bool = (type>0) .AND. (type <= ode_slv_tot_nb)
313  if (.NOT.bool) call quit(&
314  &"ode_solver_mod: ode_solver_create:&
315  & wrong argument 'type'")
316  slv%type = type
317 
318  if (present(nl_meth)) slv%NL_meth = nl_meth
319  if (present(sl_meth)) slv%SL_meth = sl_meth
320  if (present(dc_meth)) slv%DC_meth = dc_meth
321  if (present(l_meth )) slv%L_meth = l_meth
322  if (present(verb)) slv%verb = verb
323 
324  !! Opêrator splitting
325  if (present(os_meth)) then
326  slv%os = ode_opsplt(os_meth, slv%L_meth, slv%NL_meth, &
327  & check_overflow)
328  end if
329 
330  !! check overflow when solving
331  if (present(check_overflow)) then
332  slv%check_overflow = check_overflow
333  end if
334 
335  if (.NOT.valid(slv)) call quit(&
336  & "ode_solver_mod: ode_solver_create: not valid")
337 
338  end function ode_solver_create
339 
340 
343  function ode_solver_valid(slv) result(bool)
344  type(ode_solver) , intent(in) :: slv
345 
346  logical :: bool
347 
348  !! check type
349  bool = (slv%type>0) .AND. (slv%type <= ode_slv_tot_nb)
350  if (.NOT.bool) return
351 
352  !! check operator splitting
353  if (slv%type==ode_slv_os) bool = valid(slv%os)
354  if (.NOT.bool) return
355 
356  !! check ODE problem
357  bool = (slv%pb>0) .AND. (slv%pb<=ode_pb_tot_nb)
358  if (.NOT.bool) return
359 
360  !! check ODE solver depending on the ODE type
361  select case(slv%pb)
362 
363  case(ode_pb_lin)
364 
365  select case(slv%type)
366  case(ode_slv_1s)
367  bool = check_ode_method_lin_1s(slv%L_meth)
368  end select
369 
370  case(ode_pb_nl)
371 
372  select case(slv%type)
373  case(ode_slv_1s)
374 
375  bool = check_ode_method_nl_1s(slv%NL_meth)
376 
377  case(ode_slv_ms)
378 
379  bool = check_ode_method_nl_ms(slv%NL_meth)
380 
381  end select
382 
383  case(ode_pb_sl)
384 
385  select case(slv%type)
386  case(ode_slv_ms)
387 
388  bool = check_ode_method_sl_ms(slv%SL_meth)
389 
390  end select
391 
392  case(ode_pb_sl_nl)
393 
394  select case(slv%type)
395 
396  case(ode_slv_ms)
397  bool = check_ode_method_sl_ms(slv%SL_meth)
398  bool = bool .AND. check_ode_method_nl_ms(slv%NL_meth)
399 
400  case(ode_slv_dc)
401  bool = check_ode_method_sl_nl_dc(slv%DC_meth)
402  end select
403 
404  end select
405 
406  end function ode_solver_valid
407 
408 
411  subroutine ode_solver_print(slv)
412  type(ode_solver) , intent(in) :: slv
413 
414  write(*,*)"ode_solver_mod : print"
415 
416  write(*,*)" Associated ODE problem = ",&
417  & trim(name_ode_problem(slv%pb))
418 
419  write(*,*)" ODE solver type = ",&
420  & trim(name_ode_solver_type(slv%type))
421 
422  if ( valid(slv)) then
423  write(*,*)" Status = valid"
424  else
425  write(*,*)" Status = invalid"
426  end if
427 
428  select case(slv%type)
429  case(ode_slv_os)
430  call print(slv%os)
431 
432  case(ode_slv_1s)
433  if (slv%L_meth /= -1 ) then
434  write(*,*)" one-step linear solver = ",&
435  & trim(name_ode_method(slv%L_meth))
436  end if
437  if (slv%NL_meth /= -1 ) then
438  write(*,*)" one-step non-linear solver = ",&
439  & trim(name_ode_method(slv%NL_meth))
440  end if
441 
442  case(ode_slv_ms)
443  if (slv%NL_meth /= -1 ) then
444  write(*,*)" multistep non-linear solver = ",&
445  & trim(name_ode_method(slv%NL_meth))
446  end if
447  if (slv%SL_meth /= -1 ) then
448  write(*,*)" multistep semilinear solver = ",&
449  & trim(name_ode_method(slv%SL_meth))
450  end if
451 
452  case(ode_slv_dc)
453  write(*,*)" deferred correction ode solver = ",&
454  & trim(name_ode_method(slv%DC_meth))
455 
456  end select
457  write(*,*)" check overfow =", slv%check_overflow
458  write(*,*)" verbosity =", slv%verb
459 
460  end subroutine ode_solver_print
461 
462 
463 
466  function name_ode_solver_type(type) result(name)
467  integer, intent(in) :: type
468  character(len=20) :: name
469 
470  select case(type)
471  case(ode_slv_1s)
472  name="one-step"
473 
474  case(ode_slv_ms)
475  name="multistep"
476 
477  case(ode_slv_dc)
478  name="deferred corrections"
479 
480  case(ode_slv_os)
481  name="op splitting"
482 
483  case default
484  name = "invalid"
485 
486  end select
487 
488  end function name_ode_solver_type
489 
494  function check_ode_method(method, pb_type, slv_type) result(b)
495  logical :: b
496  integer, intent(in) :: method, pb_type, slv_type
497 
498  b = .false.
499  select case(pb_type)
500 
501  case(ode_pb_lin)
502 
503  select case(slv_type)
504  case(ode_slv_1s)
505  b = check_ode_method_lin_1s(method)
506 
507  end select
508 
509  case(ode_pb_nl)
510 
511  select case(slv_type)
512  case(ode_slv_ms)
513  b = check_ode_method_nl_ms(method)
514 
515  case(ode_slv_1s)
516  b = check_ode_method_nl_1s(method)
517 
518  end select
519 
520  case(ode_pb_sl)
521 
522  select case(slv_type)
523  case(ode_slv_ms)
524  b = check_ode_method_sl_ms(method)
525 
526  end select
527 
528  case(ode_pb_sl_nl)
529 
530  select case(slv_type)
531  case(ode_slv_dc)
532  b = check_ode_method_sl_nl_dc(method)
533 
534  end select
535 
536  end select
537 
538  end function check_ode_method
539 
540 
541 
552  function ode_solution_create(slv, pb) result(sol)
553  type(ode_solution) :: sol
554  type(ode_solver) , intent(in) :: slv
555  type(ode_problem) , intent(in) :: pb
556 
557  logical :: bool
558 
559  if (choral_verb>0) write(*,*) &
560  &"ode_solver_mod : ode_solution_create"
561 
562  bool = valid(pb)
563  if (.NOT.bool) call quit(&
564  & "ode_solver_mod: ode_solution_create: 'pb' not valid")
565 
566  bool = valid(slv)
567  if (.NOT.bool) call quit(&
568  & "ode_solver_mod: ode_solution_create: 'slv' not valid")
569 
570  !! define ODE solution
571  select case(pb%type)
572 
573  case(ode_pb_lin)
574 
575  select case(slv%type)
576  case(ode_slv_1s)
577  call create_ode_lin_1s_sol(sol, pb, slv%L_meth)
578 
579  end select
580 
581  case(ode_pb_nl)
582 
583  select case(slv%type)
584  case(ode_slv_1s)
585  call create_ode_nl_1s_sol(sol, pb, slv%NL_meth)
586 
587  case(ode_slv_ms)
588  call create_ode_nl_ms_sol(sol, pb, slv%NL_meth)
589 
590  end select
591 
592  case(ode_pb_sl)
593 
594  select case(slv%type)
595  case(ode_slv_os)
596  call create_ode_opsplt_sol(sol, pb, slv%os)
597 
598  case(ode_slv_ms)
599  call create_ode_sl_ms_sol(sol, pb, slv%SL_meth)
600 
601  end select
602 
603  case(ode_pb_sl_nl)
604 
605  select case(slv%type)
606  case(ode_slv_os)
607  call create_ode_opsplt_sol(sol, pb, slv%os)
608 
609  case(ode_slv_ms)
610  call create_ode_sl_nl_ms_sol(sol, pb, &
611  & slv%SL_meth, slv%NL_meth)
612 
613  case(ode_slv_dc)
614  call create_ode_sl_nl_dc_sol(sol, pb, slv%DC_meth)
615 
616  end select
617 
618  end select
619 
620  end function ode_solution_create
621 
622 
625  subroutine homogeneous_initialcond(sol, pb, slv, t0, y0)
626  type(ode_solution) , intent(inout):: sol
627  type(ode_problem ) , intent(in) :: pb
628  type(ode_solver) , intent(in) :: slv
629  real(RP) , intent(in) :: t0
630  real(RP), dimension(:), intent(in) :: y0
631 
632  integer :: jj, ii, N, Na
633 
634  if (choral_verb>2) write(*,*) &
635  &"ode_solver_mod :&
636  & initialCond = homogeneous"
637 
638  if (size(y0)/=sol%N) call quit( &
639  & "ode_solver_mod: homogeneous_initialCond:&
640  & wrong size for Y0" )
641 
642  n = sol%N
643  na = sol%Na
644 
645  associate( y=>sol%Y, by=>sol%BY, ay=>sol%AY, &
646  & x=>pb%X, v=>sol%V)
647 
648  !! initialise Y
649  do ii=1, sol%dof
650  do jj=1, sol%nY
651  y(:,jj,ii) = y0(:)
652  end do
653  end do
654 
655  !! initialise V
656  do jj=1, sol%nV
657  v(:,jj) = y0(n)
658  end do
659 
660  !! initialise AY and BY
661  if (sol%NFY==0) return
662 
663  do ii=1, sol%dof
664  call pb%AB(ay(:,1,ii), by(:,1,ii), x(:,ii), &
665  & t0, y0, n, na)
666  end do
667 
668  do jj=2, sol%nFY
669  ay(:,jj,:) = ay(:,1,:)
670  by(:,jj,:) = by(:,1,:)
671  end do
672 
673  if (na == n) then
674  if ( (pb%type==ode_pb_sl).OR.(pb%type==ode_pb_sl_nl) ) then
675  if (slv%type==ode_slv_ms) then
676 
677  do ii=1, sol%dof
678  do jj=1, sol%nFY
679 
680  by(n,jj,ii) = by(n,jj,ii) &
681  & + ay(n,jj,ii)*y0(n)
682  end do
683  end do
684 
685  end if
686  end if
687  end if
688 
689  end associate
690 
691  end subroutine homogeneous_initialcond
692 
693 
714  subroutine ode_solver_solve(sol, slv, pb, t0, T, dt, KInv, output)
715  type(ode_solution) , intent(inout) :: sol
716  type(ode_problem) , intent(in) :: pb
717  type(ode_solver) , intent(inout) :: slv
718  real(RP) , intent(in) :: dt, t0, T
719 
720  procedure(linsystem_solver), optional :: KInv
721  type(ode_output) , optional :: output
722 
723  procedure(ode_output_proc), pointer :: out => null()
724  logical :: bool
725  real(RP) :: time
726  real(DP) :: cpu
727 
728  if (slv%verb>0) write(*,*) &
729  &"ode_solver_mod : solve dt =", dt
730 
731  bool = valid(pb)
732  if (.NOT.bool) call quit(&
733  &"ode_solver_mod: ode_solver_solve: arg 'pb' not valid")
734 
735  bool = valid(slv)
736  if (.NOT.bool) call quit(&
737  &"ode_solver_mod: ode_solver_solve: arg 'slv' not valid")
738 
739  !! output defined with 'out'
740  !!
741  if (present(output)) then
742 
743  !! check
744  if( .NOT.output%assembled &
745  & .OR. (output%nbDof /= pb%dof) &
746  & .OR. (output%N /= pb%N ) ) then
747  call quit(&
748  & "ode_solver_mod: ode_solver_solve: arg 'out' not valid")
749  end if
750 
751  out => local_output
752 
753  else
754  out => slv%output
755 
756  end if
757 
758 
759 
760  cpu = clock()
761  sol%ierr = 0
762 
763  select case(pb%type)
764 
765  case(ode_pb_lin)
766  select case(slv%type)
767  case(ode_slv_1s)
768 
769  call solve_ode_lin_1s(sol, pb, t0, t, dt, &
770  & slv%L_meth, out, kinv, slv%kry)
771 
772  end select
773 
774  case(ode_pb_nl)
775  select case(slv%type)
776  case(ode_slv_1s)
777 
778  call solve_ode_nl_1s(sol, pb, t0, t, dt, &
779  & slv%NL_meth, out, &
780  & slv%check_overflow)
781 
782  case(ode_slv_ms)
783 
784  call solve_ode_nl_ms(sol, pb, t0, t, dt, &
785  & slv%NL_meth, out, &
786  & slv%check_overflow)
787 
788  end select
789 
790  case(ode_pb_sl)
791  select case(slv%type)
792  case(ode_slv_ms)
793 
794  call solve_ode_sl_ms(sol, pb, t0, t, dt, &
795  & slv%SL_meth, out, &
796  & slv%check_overflow, kinv, slv%kry)
797 
798  case(ode_slv_os)
799  call solve_ode_opsplt(sol, t0, t, dt, &
800  & slv%os, pb, slv%kry, out)
801  end select
802 
803  case(ode_pb_sl_nl)
804 
805  select case(slv%type)
806  case(ode_slv_ms)
807 
808  call solve_ode_sl_nl_ms(sol, pb, t0, t, dt, &
809  & slv%SL_meth, slv%NL_meth, out, &
810  & slv%check_overflow, kinv, slv%kry)
811 
812  case(ode_slv_os)
813  call solve_ode_opsplt(sol, t0, t, dt, &
814  &slv%os, pb, slv%kry, out)
815 
816  case(ode_slv_dc)
817 
818  call solve_ode_sl_nl_dc(sol, pb, t0, t, dt, &
819  & slv%DC_meth, out, &
820  & slv%check_overflow, kinv, slv%kry)
821  end select
822 
823  end select
824 
825  if (slv%verb>1) then
826  cpu = clock() - cpu
827  write(*,*)"ode_solver_mod : solve end, CPU =",&
828  & real(cpu, sp)
829  end if
830 
831  if (sol%ierr/=0) then
832  if (slv%verb>0) write(*,*) &
833  &"ode_solver_mod :&
834  & solve = ERROR DETECTED, IERR = ", sol%ierr
835 
836  else
837  !! external post treatment
838  !! special output at computation end
839  time = t
840  bool = .true.
841  call out(time, sol, bool)
842  end if
843 
844  contains
845 
846  subroutine local_output(tn, s, stop)
847  real(RP) , intent(in) :: tn
848  type(ode_solution), intent(in) :: s
849  logical , intent(inout) :: stop
850 
851  call output_proc(output, tn, s, stop)
852 
853  end subroutine local_output
854 
855  end subroutine ode_solver_solve
856 
857 
858 
871  subroutine set_ode_solver_output(slv, output)
872  type(ode_solver), intent(inout) :: slv
873  procedure(ode_output_proc) :: output
874 
875  slv%output => output
876 
877  end subroutine set_ode_solver_output
878 
879 end module ode_solver_mod
DEFERRED CORRECTION SOLVERS FOR A SEMILINEAR ODE coupled with a non-linear ODE system ...
MULTISTEP SOLVERS FOR A SEMILINEAR ODE coupled with a non-linear ODE system
type(ode_solver) function ode_solver_create(pb, type, os_meth, SL_meth, NL_meth, L_meth, DC_meth, check_overflow, verb)
Constructor for the type ode_solver
subroutine homogeneous_initialcond(sol, pb, slv, t0, y0)
homogeneous initial condition
BASIC TOOLS
Definition: basic_tools.F90:8
subroutine ode_solver_clear(slv)
destructor for ode_solver
double precision function, public clock()
Returns the internal clock time
Definition: basic_tools.F90:94
DERIVED TYPE ode_problem: definition of ODE/PDE problems
subroutine, public create_ode_sl_nl_ms_sol(sol, pb, SL_meth, NL_meth)
Create the solution data structure.
subroutine, public solve_ode_opsplt(sol, t0, T, dt, os, pb, kry, out)
solve : operator splitting with constant time step
subroutine, public create_ode_sl_nl_dc_sol(sol, pb, method)
Create the solution data structure.
subroutine, public create_ode_opsplt_sol(sol, pb, os)
The type krylov defines the settings of a linear solver.
Definition: krylov_mod.f90:63
subroutine, public solve_ode_sl_nl_dc(sol, pb, t0, T, dt, meth, out, check_overflow, Kinv, kry)
solve : multistep with constant time step
subroutine, public quit(message)
Stop program execution, display an error messahe.
subroutine, public create_ode_sl_ms_sol(sol, pb, method)
create memory for the ode_solution structure &#39;sol&#39;
integer, parameter, public sp
simple precision for real numbers
Definition: real_type.F90:60
integer, parameter ode_pb_tot_nb
Number of ode problem types.
logical function, public check_ode_method_nl_ms(method)
is &#39;method&#39; a multi-step non-linear ODE solver ?
The type opSplt defines operator spltting methods.
integer, parameter ode_slv_ms
multistep
logical function, public check_ode_method_sl_nl_dc(method)
is &#39;method&#39; a DC method for SL_NL problem ?
character(len=15) function, public name_ode_method(method)
Get ODE method name.
Definition: ode_def.f90:68
subroutine, public solve_ode_sl_nl_ms(sol, pb, t0, T, dt, SL_meth, NL_meth, out, check_overflow, Kinv, kry)
solve : multistep with constant time step
REAL NUMBERS PRECISION IN CHORAL: selects simple/double/quad
Definition: real_type.F90:33
subroutine, public solve_ode_lin_1s(sol, pb, t0, T, dt, method, out, Kinv, kry)
solve with constant time-step
logical function, public check_ode_method_nl_1s(method)
is &#39;method&#39; a one-step non-linear ODE solver ?
integer, parameter ode_pb_sl
SemiLinear ODE : with .
ONE-STEP SOLVERS FOR LINEAR ODEs
integer, parameter ode_pb_lin
Linear ODE : .
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.
real(rp) function e(x, v1, v2)
CHORAL CONSTANTS
logical function, public check_ode_method(method, pb_type, slv_type)
check whether the ode method &#39;method&#39; is available for the problem type &#39;pb_type&#39; and for the solver ...
subroutine, public solve_ode_nl_1s(sol, pb, t0, T, dt, method, out, check_overflow)
solve with constant time-step
integer, parameter ode_slv_dc
deferred corrections
subroutine, public output_proc(out, time, s, final)
Pre-defined output for ODE resolution
integer, parameter ode_slv_tot_nb
number of ode solver types
DERIVED TYPE ode_output: handles output for PDE/ODE simulations
logical function, public check_ode_method_sl_ms(method)
is &#39;method&#39; a multi-step semilinear ODE solver ?
subroutine, public solve_ode_nl_ms(sol, pb, t0, T, dt, method, out, check_overflow)
solve with constant time-step
set initial conditions
subroutine ode_solver_solve(sol, slv, pb, t0, T, dt, KInv, output)
Solve an ODE with constant time step
BOTTOM LEVEL MODULE FOR ODEs
Definition: ode_def.f90:12
TOP LEVEL MODULE FOR ODEs, derived type ode_solver
integer choral_verb
Verbosity level.
subroutine, public create_ode_nl_1s_sol(sol, pb, method)
create memory for the ode_solution structure &#39;sol&#39;
integer, parameter ode_pb_sl_nl
SemiLinear ODE coupled with a non-linear ODE system for with .
character(len=20) function, public name_ode_solver_type(type)
name the type of ode_solver
integer, parameter ode_slv_os
operator splitting
The type ode_output handles output for ODE simulations.
subroutine local_output(tn, s, stop)
ONE-STEP SOLVERS FOR NON LINEAR ODEs
subroutine, public create_ode_lin_1s_sol(sol, pb, method)
allocate memory for the ode_solution structure &#39;sol&#39;
subroutine, public set_ode_solver_output(slv, output)
Load a user defined output for ODE resolution
MULTISTEP SOLVERS FOR SEMILINEAR ODEs
DEFINITION OF GLOBAL VARIABLES FOR THE LIBRARY CHORAL
subroutine, public solve_ode_sl_ms(sol, pb, t0, T, dt, method, out, check_overflow, Kinv, kry)
solve with constant time-step Case where pbN = 1
logical function, public check_ode_method_lin_1s(method)
is &#39;method&#39; a one-step linear ODE solver ?
character(len=15) function, public name_ode_problem(type)
ode_problem name
MULTISTEP SOLVERS FOR NON LINEAR ODEs
DERIVED TYPE ode_opSplt: operator splitting methods for ODEs.
Generic &#39;solve&#39; for ODEs.
Type ode_problem: definition of ODE/PDE problems.
integer, parameter ode_pb_nl
Non-Linear ODE system: , .
subroutine ode_solver_print(slv)
print ode_solver
DERIVED TYPE krylov: for the resolution of linear systems
Definition: krylov_mod.f90:25
subroutine, public void_ode_output(tn, s, stop)
void output for ode resolution
integer, parameter ode_slv_1s
onestep
integer, parameter kry_cg
CG linear solver.
logical function ode_solver_valid(slv)
check ode_solver
subroutine, public create_ode_nl_ms_sol(sol, pb, method)
create memory for the ode_solution structure &#39;sol&#39;